A three-dimensional image post-processing method for CBCT beam hardening correction

Abstract

This invention provides a three-dimensional image post-processing method for CBCT beam hardening correction, comprising: performing three-dimensional voxel density layering on CBCT three-dimensional reconstructed volume data to generate a three-dimensional density layering matrix, and calculating the beam path cumulative density of each voxel; calculating adaptive correction coefficients on the CBCT three-dimensional reconstructed data based on the three-dimensional density layering matrix and the path cumulative density, and performing CT value adaptive correction to obtain preliminary correction data; and performing artifact residue detection and secondary optimization on the preliminary correction data for the tumor target area layer and low-contrast soft tissue layer to obtain optimized correction data. The method of this invention solves the technical problems of existing CBCT beam hardening correction techniques, such as poor artifact removal, CT value distortion, neglect of three-dimensional correlation, inability to adapt to real-time radiotherapy workflows, easy loss of target area and endangered organ details, and difficulty in clinical integration.

Description

Technical Field

[0001] This invention belongs to the field of medical image processing technology, specifically relating to a three-dimensional image post-processing method for CBCT beam hardening correction, which can effectively eliminate cup artifacts, streak artifacts and CT value distortion, accurately restore the imaging quality of tumor target area and organs at risk, and improve the accuracy of radiotherapy positioning, target delineation and dose calculation. Background Technology

[0002] CBCT, due to its low radiation, fast imaging speed, and ability to perform real-time bedside scanning, has become a core device for image-guided localization in tumor radiotherapy. However, X-ray beam hardening can cause problems in CBCT, such as grayscale distortion around high-density tissues, streak artifacts around metal implants, and target area CT value shifts. These issues directly lead to target area delineation deviations, blurred boundaries of organs at risk, positioning and registration errors, and inaccurate dose calculations, severely affecting the precision of radiotherapy. Conventional 3D filtering only achieves noise smoothing and lacks a targeted beam hardening correction mechanism, easily resulting in the loss of tumor target area and details of organs at risk.

[0003] Existing calibration techniques have significant drawbacks:

[0004] 1. Traditional water model / polynomial correction is a global static correction, which cannot adapt to the complex tissue density heterogeneity of tumor patients and has poor artifact removal effect;

[0005] 2. Deep learning correction relies on a large amount of labeled data, has high computing power requirements, and weak generalization ability, making it difficult to integrate into bedside CBCT systems for radiotherapy linear accelerators and unable to meet real-time positioning requirements.

[0006] Therefore, there is an urgent need for a beam hardening correction method that is compatible with CBCT three-dimensional volume data, computationally efficient, easy to integrate, and specifically designed for tumor radiotherapy localization, in order to solve the technical problems of artifact removal, CT value fidelity, preservation of target / organs at risk details, and clinical real-time compatibility. Summary of the Invention

[0007] The purpose of this invention is to provide a three-dimensional image post-processing method for CBCT beam hardening correction, in order to solve the technical problems of existing CBCT beam hardening correction technology, such as poor artifact removal, CT value distortion, neglect of three-dimensional correlation, inability to adapt to real-time radiotherapy workflow, easy loss of target area and endangered organ details, and difficulty in clinical integration.

[0008] To achieve the above objectives, the present invention provides a three-dimensional image post-processing method for CBCT beam hardening correction, comprising:

[0009] The CBCT 3D reconstructed volume data is subjected to 3D voxel density layering to generate a 3D density layering matrix, and the beam path cumulative density of each voxel is calculated.

[0010] Based on the three-dimensional density layer matrix and path cumulative density, adaptive correction coefficients are calculated for CBCT three-dimensional reconstruction data, and adaptive correction of CT values ​​is performed to obtain preliminary correction data.

[0011] The preliminary correction data were subjected to artifact residue detection and secondary optimization in the tumor target layer and low-contrast soft tissue layer to obtain optimized correction data.

[0012] The three-dimensional density layer matrix DensityLayer includes a density layer label for each voxel, which is determined based on the CT value physical characteristics of key tissues in tumor radiotherapy and the clinical needs of tumor radiotherapy localization.

[0013] Based on the physical characteristics of CT values ​​of key tissues in tumor radiotherapy, three-dimensional voxels are divided into five density levels, including air / cavity layer, low-contrast soft tissue layer, tumor target area layer, bone layer, and metal layer. The tumor target area layer is used to meet the precise positioning requirements of tumor radiotherapy.

[0014] The CT values ​​for each density level are as follows: Air / cavity layer: CT < 500 HU; Low-contrast soft tissue layer: -500 HU < CT ≤ 200 HU; Tumor target layer: 200 HU < CT ≤ 300 HU; Bone layer: 300 HU < CT ≤ 2000 HU; Metal layer: 2000 HU < CT.

[0015] The path cumulative density of each voxel is calculated to form a three-dimensional density layer matrix. Specifically, this includes: calculating the X-ray incident path of each voxel P(xx,yy,zz) based on the X-ray source position, detector parameters, and patient positioning in the CBCT scan, and calculating the path cumulative density D. sum (xx,yy,zz);

[0016] Path cumulative density D sum The formula for calculating (xx,yy,zz) is:

[0017] ,

[0018] in, Let the coordinates be those of a voxel. The total number of voxels along the X-ray incident path, CT i The first X-ray in the incident path Raw CT values ​​of individual pixels, w i For the first The density hierarchy weighting coefficient of a voxel, where i is the ordinal number of the voxel.

[0019] The larger the adaptive correction coefficient, the smaller the correction degree; the adaptive correction coefficient is related to the density level, and the adaptive correction coefficient corresponding to the density level of the tumor target area layer and the low contrast soft tissue layer is higher than the adaptive correction coefficient of other density levels; and the adaptive correction coefficient is negatively correlated with the path cumulative density.

[0020] The adaptive correction coefficients k(xx,yy,zz) are:

[0021] ,

[0022] Where k0(L) is the baseline correction coefficient for the density level, and L is the density level label. D is the attenuation adjustment coefficient. max This represents the maximum value of the path cumulative density. Let the coordinates be those of a voxel. For path cumulative density;

[0023] Preliminary correction data for:

[0024] ,

[0025] in, These are the original CT values. For preliminary data correction, β is the adaptive correction coefficient, and β is the CT value offset compensation amount.

[0026] When the density level label is air / cavity layer, k0(L) is When the density level label is a low-contrast soft tissue layer, k0(L) is... When the density level label is the tumor target layer, k0(L) is... When the density level label is the bone layer, k0(L) is... When the density level label is a metal layer, k0(L) is... .

[0027] The preliminary corrected data underwent artifact residue detection and secondary optimization in the tumor target layer and low-contrast soft tissue layer to obtain optimized corrected data, specifically including:

[0028] Based on the preliminary correction data, calculate the three-dimensional gradient value and detect and identify artifact regions based on the three-dimensional gradient value;

[0029] For the artifact residue regions marked in the artifact mask matrix, fine correction is performed using three-dimensional neighborhood mean smoothing and strong constraints of radiotherapy-specific density levels to obtain optimized correction data.

[0030] Based on the preliminary correction data, the three-dimensional gradient value is calculated, and artifact regions are detected and identified based on the three-dimensional gradient value. Specifically, this includes:

[0031] A 3D neighborhood window is set for each voxel P(xx,yy,zz), and an artifact detection gradient threshold G is pre-set. th ;

[0032] Based on the preliminary correction data and the three-dimensional neighborhood window, the three-dimensional gradient value G(xx,yy,zz) of the voxel is calculated to characterize the gray-level change rate between the voxel and its neighborhood.

[0033] The three-dimensional gradient value G(xx,yy,zz) is compared with the artifact detection gradient threshold G. th The comparison is performed, and the artifact residue area is determined based on the comparison results; the artifact residue areas around the tumor target layer and the low contrast soft tissue layer in the three-dimensional density layering matrix are marked to form an artifact mask matrix;

[0034] Optimize calibration data for:

[0035] ,

[0036] Where N is the 3×3×3 effective neighborhood of the voxel being corrected, and M is the number of voxels in the effective neighborhood. It is a three-dimensional artifact mask matrix. For the initial data correction, γ(L) is the smoothing constraint coefficient of the density level, Q is the voxel, L is the density level label, and (xx,yy,zz) are the voxel coordinates.

[0037] After obtaining the optimized correction data, the process also includes outputting the optimized correction data, specifically: converting the optimized correction data into a radiotherapy-specific format and outputting it in a format suitable for clinical use in tumor radiotherapy.

[0038] The adaptive beam hardening correction method based on three-dimensional voxel density stratification of this invention adopts a three-level three-dimensional cascaded process of "three-dimensional density stratification and beam path analysis → radiotherapy-specific adaptive stratification correction → secondary optimization of target area / organ at risk artifacts". It removes artifacts and restores CT values ​​while preserving tumor and organ at risk details. It can be directly integrated into radiotherapy CBCT systems to meet the needs of precise tumor radiotherapy. It solves the technical problems of existing CBCT beam hardening correction technology, such as poor artifact removal effect, CT value distortion, neglect of three-dimensional correlation, inability to adapt to real-time radiotherapy process, easy loss of target area and organ at risk details, and high difficulty in clinical integration. Attached Figure Description

[0039] Figure 1 This is a flowchart of a three-dimensional image post-processing method for CBCT beam hardening correction according to the present invention.

[0040] Figure 2AImage of a uniform water layer on a Catphan700 phantom before correction;

[0041] Figure 2B The image of a uniform water layer from a Catphan700 phantom test image is corrected using the method of this invention.

[0042] Figure 3A This is an image of a slice layer of a head phantom after correction using the method of the present invention;

[0043] Figure 3B yes Figure 3A The image before correction;

[0044] Figure 3C This is an image of another head phantom slice layer after correction using the method of the present invention;

[0045] Figure 3D yes Figure 3C The image before correction. Detailed Implementation

[0046] The following specific examples illustrate the implementation of the present invention. Those skilled in the art can easily understand other advantages and effects of the present invention from the content disclosed in this specification. The present invention can also be implemented or applied through other different specific embodiments, and various details in this specification can also be modified or changed based on different viewpoints and applications without departing from the spirit of the present invention. It should be noted that, unless otherwise specified, the following embodiments and features described therein can be combined with each other.

[0047] It should be noted that the illustrations provided in the following embodiments are only schematic representations of the basic concept of the present invention. Therefore, the drawings only show the components related to the present invention and are not drawn according to the actual number, shape and size of the components in the actual implementation. In the actual implementation, the form, quantity and proportion of each component can be arbitrarily changed, and the layout of the components may also be more complex.

[0048] The principle of the three-dimensional image post-processing method for CBCT beam hardening correction of the present invention is as follows: for the three-dimensional volume data of CBCT for tumor treatment, radiotherapy-specific layered adaptive correction and precise optimization of target area artifacts are performed. Based on the density characteristics of three-dimensional voxels and the spatial correlation of X-ray beam paths, combined with the requirements of precise identification of target area and organs at risk for tumor radiotherapy positioning, a full-process three-dimensional operation is designed to achieve a high degree of fit between beam hardening correction and tumor radiotherapy positioning needs.

[0049] Taking image-guided localization using a head tumor radiotherapy phantom as an example, the three-dimensional image post-processing method for CBCT beam hardening correction of the present invention specifically includes:

[0050] Step S1: Acquire CBCT three-dimensional reconstruction data for tumor treatment;

[0051] The process involves collecting raw CBCT scan projection data from the radiotherapy linear accelerator of tumor radiotherapy patients, and then using a radiotherapy-specific FDK reconstruction algorithm to obtain CBCT 3D reconstruction data. The CBCT 3D reconstruction data is in the form of a 3D data matrix InputData[X×Y×Z], where X, Y, and Z represent the pixel dimensions of the CBCT 3D reconstruction data in the transverse, coronal, and sagittal planes, respectively. The data format of the 3D data matrix is ​​floating-point, and the voxel values ​​are the original CT values ​​(HU).

[0052] The CBCT 3D reconstruction data is directly interfaced with the image guidance system of the radiotherapy linear accelerator to prepare for subsequent 3D beam hardening correction.

[0053] In this embodiment, the head radiotherapy simulation phantom is scanned by bedside CBCT using a radiotherapy linear accelerator to acquire raw projection data. The raw data is then reconstructed using a radiotherapy-specific FDK algorithm to obtain 512×512×512-dimensional floating-point three-dimensional volume data. The original CT values ​​of the voxels range from -1000 to 3000 HU (including nasal air, brain tissue, simulated tumor target area, skull base bones, and metal positioning markers). This data is used as the input data matrix InputData and directly connected to the accelerator image guidance system.

[0054] Step S2: Perform density layering on the CBCT 3D reconstruction data to generate a 3D density layering matrix DensityLayer, and calculate the path cumulative density of each voxel;

[0055] The three-dimensional density layer matrix DensityLayer includes a density layer label for each voxel, which is determined based on the physical characteristics of CT values ​​of key tissues in tumor radiotherapy and the clinical needs of tumor radiotherapy localization.

[0056] In this embodiment, based on the physical characteristics of CT values ​​of key tissues in tumor radiotherapy, the three-dimensional voxels are divided into five density levels: air / cavity layer, low-contrast soft tissue layer, tumor target area layer, bone layer, and metal layer. The newly added tumor target area layer is used to adapt to the precise positioning requirements of tumor radiotherapy. The CT value threshold for each density level can be flexibly adjusted according to different tumor locations such as lung cancer, head and neck tumors, and pelvic tumors. The CT values ​​for each density level are as follows:

[0057] Air / cavitary layer: CT < 500 HU (e.g., intrapulmonary cavities, intestinal gas), density level label is 1;

[0058] Low-contrast soft tissue layers: -500 HU < CT ≤ 200 HU (such as normal soft tissues like lung, liver, and pancreas), density level label is 2;

[0059] Tumor target layer: 200 HU < CT ≤ 300 HU (such as most solid tumor lesions, including small low-contrast lesions), density level label is 3;

[0060] Bone layer: 300 HU < CT ≤ 2000 HU (such as bony landmarks like vertebrae, sternum, and pelvis), density level label is 4;

[0061] Metallic layer: 2000 HU < CT (e.g., radiotherapy metal positioning markers, in vivo metal implants), density level label is 5.

[0062] For each voxel P(xx,yy,zz), a corresponding density layer label L(xx,yy,zz)∈{1,2,3,4,5} is assigned to form a three-dimensional density layer matrix DensityLayer, whose data form is [X×Y×Z], to achieve precise stratification of the tumor target area and surrounding tissues.

[0063] The path cumulative density of each voxel is calculated to form a three-dimensional density layer matrix. Specifically, this includes: calculating the X-ray incident path of each voxel P(xx,yy,zz) based on the X-ray source position, detector parameters, and patient positioning in the CBCT scan, and calculating the path cumulative density D. sum (xx,yy,zz), thus providing a basis for subsequent adaptive layered correction. The X-ray incident path is usually a straight line.

[0064] The X-ray incident path is set to a step size of 1 voxel (to meet the high-precision requirements of radiotherapy positioning). All voxels along the X-ray incident path are traversed, and the density is calculated by summing the values ​​based on a density of 1 voxel unit. The cumulative density D along the path is then calculated. sum The formula for calculating (xx,yy,zz) is:

[0065] ,

[0066] in, Let the coordinates be those of a voxel. The total number of voxels along the X-ray incident path, CT i The first X-ray in the incident path Raw CT values ​​of individual pixels, w i For the first The density hierarchy weighting coefficient of individual elements is set to 1.0 for the tumor target layer to preserve its density characteristics to the maximum extent. The weighting coefficient can be adjusted according to different tumor locations.

[0067] From the above formula, we can see that the path cumulative density D of all voxels P(xx,yy,zz) located on the same X-ray incident path is... sum (xx,yy,zz) should be the same.

[0068] In this embodiment, the first When the density hierarchy label of an individual element is air / cavity layer, the first Individual density hierarchy weighting coefficient w i = 0.1; When it is a low-contrast soft tissue layer, the first Individual density hierarchy weighting coefficient w i = 0.8; when it is the tumor target layer, the first Individual density hierarchy weighting coefficient w i = 1.0; when it is the bone layer, the first Individual density hierarchy weighting coefficient w i = 1.5; when it is a metal layer, the first Individual density hierarchy weighting coefficient w i = 3.0.

[0069] In this embodiment, density layers are defined according to radiotherapy-specific thresholds (air / cavity < -500 HU, soft tissues such as brain tissue -500~200 HU, simulated tumor target area 200~300 HU, skull base / cervical spine 300~2000 HU, metal positioning markers >2000 HU), generating a three-dimensional density layer matrix DensityLayer to accurately mark the simulated tumor target area and the brainstem and spinal cord (organs at risk); the X-ray incident path step size is set to 1 voxel, and the path cumulative density D of each voxel is calculated. sum (xx,yy,zz), where the path cumulative density D sum The maximum value of (xx,yy,zz) is 2600 HU.

[0070] Step S3: Based on the three-dimensional density layer matrix and path cumulative density D sum Adaptive correction coefficients are calculated for CBCT three-dimensional reconstruction data, and adaptive correction of CT values ​​is performed to obtain preliminary correction data;

[0071] The adaptive correction coefficients are derived based on the exponential law of X-ray attenuation. Furthermore, the correction degree for CT values ​​in the density layers of the tumor target area and low-contrast soft tissue layers is less than that in other density layers. This results in greater protection of the CT value accuracy of the tumor target area and the layers containing the main organs at risk compared to other regions, avoiding over-correction that leads to loss of detail and suppressing beam hardening artifacts. Thus, a three-dimensional adaptive beam hardening layered correction specifically for tumor radiotherapy is achieved, enabling the correction coefficients to dynamically adapt to the degree of beam hardening while providing correction and protection for the tumor target area and organs at risk.

[0072] In this embodiment, a larger adaptive correction coefficient results in a smaller correction degree. The adaptive correction coefficient is related to the density level; the adaptive correction coefficients corresponding to the density levels of the tumor target area layer and the low-contrast soft tissue layer (the layer where the organs at risk are mainly located) are higher than those of other density levels. Furthermore, the adaptive correction coefficient is negatively correlated with the path cumulative density.

[0073] The adaptive correction coefficient k(xx,yy,zz) for each voxel is:

[0074] ,

[0075] Where k0(L) is the baseline correction coefficient for the density level, and L is the density level label. D is the attenuation adjustment coefficient. max This represents the maximum value of the path cumulative density. Let the coordinates be those of a voxel. This represents the cumulative density of the path.

[0076] In this embodiment, the reference correction coefficient k0(L) for the density level is: when the density level label is an air / cavity layer, k0(L) is... When the density level label is a low-contrast soft tissue layer, k0(L) is... When the density level label is the tumor target layer, k0(L) is... The reason why the correction coefficient for the tumor target layer is slightly larger than that for soft tissue is that tumor cells are more densely packed, have abundant stroma, and different water content than normal muscle / fat. Therefore, the CT value is usually higher than that of normal soft tissue, and the corresponding correction coefficient will be slightly higher than that of normal soft tissue areas. When the density level label is bone layer, k0(L) is... When the density level label is a metal layer, k0(L) is... Attenuation adjustment coefficient Default is To reduce the decay rate and protect target details.

[0077] Specifically, based on the adaptive correction coefficient, the original CT value of each voxel in the CBCT three-dimensional reconstruction data is adaptively corrected to obtain preliminary correction data.

[0078] In this embodiment, preliminary correction data For three-dimensional volume data, the calculation formula is:

[0079] ,

[0080] Among them, CT raw The original CT values, CT pre For preliminary data correction, β is the adaptive correction coefficient, and β is the CT value offset compensation amount. The default CT value offset compensation amount is 5HU, which is used to slightly compensate for the overall grayscale offset and restore the actual density characteristics of the tumor target area and the surrounding normal tissue to the greatest extent.

[0081] In this embodiment, the baseline correction coefficient k0 for each radiotherapy level, the attenuation adjustment coefficient α=0.5, and the CT value offset compensation amount β=5 HU are set. The adaptive correction coefficient k (xx,yy,zz) for each voxel is calculated to correct the original CT value and obtain preliminary correction data. At this time, the cup-shaped artifact next to the skull base bones and the stripe artifact next to the metal positioning marks have been initially suppressed, and the simulated tumor target area CT value has no obvious distortion.

[0082] Step S4: Perform artifact residue detection and secondary optimization on the preliminary correction data for the tumor target area layer and low-contrast soft tissue layer to obtain optimized correction data; thereby removing artifact residues around the tumor target area and organs at risk that may still exist after preliminary correction and have the greatest impact on radiotherapy localization.

[0083] Step S4 specifically includes:

[0084] Step S41: Calculate the three-dimensional gradient value based on the preliminary correction data and detect and identify artifact regions based on the three-dimensional gradient value;

[0085] Step S41 specifically includes:

[0086] Step S411: Set a 3D neighborhood window for each voxel P(xx,yy,zz), and pre-set the artifact detection gradient threshold G. th ;

[0087] In this embodiment, the size of the 3D neighborhood window is 3×3×3. The artifact detection gradient threshold G... th The default value is 60 HU / voxel, which is lower than the normal threshold, to avoid missing slight artifacts around the target area.

[0088] Step S412: Based on the preliminary correction data and the three-dimensional neighborhood window, calculate the three-dimensional gradient value G(xx,yy,zz) of the voxel to characterize the gray-level change rate between the voxel and its neighborhood.

[0089] The formula for calculating the three-dimensional gradient value G(xx,yy,zz) is:

[0090] ,

[0091] in, For preliminary data correction, Let x, y, and z be the coordinates of the voxel, and x, y, and z be the coordinates in three-dimensional space.

[0092] Step S413: Combine the 3D gradient value G(xx,yy,zz) with the artifact detection gradient threshold G th The comparison is performed, and the artifact residue area is determined based on the comparison results. The artifact residue areas around the tumor target layer and the low-contrast soft tissue layer in the three-dimensional density layer matrix DensityLayer are marked to form an artifact mask matrix.

[0093] The artifact mask matrix is ​​a three-dimensional artifact mask matrix ArtifactMask[X×Y×Z]. The value of the artifact region is 1, and the value of the normal region is 0.

[0094] If G(xx,yy,zz) > G th If the region is found to be a residual artifact region, it is determined to be such a region. The artifact detection gradient threshold G... th The default value is 60 HU / voxel, which is lower than the usual threshold to avoid missing slight artifacts around the target area.

[0095] The artifact residue region around the tumor target layer and the low-contrast soft tissue layer refers to the 3x3x3 neighborhood of voxels that satisfy the density level label L(xx,yy,zz)∈{2,3} in the three-dimensional density layering matrix DensityLayer.

[0096] Step S42: For the artifact residue regions marked in the artifact mask matrix, fine correction is performed using three-dimensional neighborhood mean smoothing and strong constraints of radiotherapy-specific density levels to obtain optimized correction data. Since fine correction is performed only in the artifact residue regions around the tumor target layer and the low-contrast soft tissue layer, it can strictly avoid the loss of edge details between the tumor target layer and the low-contrast soft tissue layer.

[0097] Optimize calibration data for:

[0098] ,

[0099] Where N is the 3×3×3 effective neighborhood of the voxel being corrected, and M is the number of voxels in the effective neighborhood. It is a three-dimensional artifact mask matrix. For initial data correction, γ(L) is the smoothing constraint coefficient of the density level, Q is the voxel, L is the density level label, and (xx,yy,zz) are the voxel coordinates. The larger γ(L) is, the stronger the smoothing force. The smoothing constraint coefficient of the density level of the tumor target layer is set to the smallest value except for the metal layer, thereby minimizing the impact of the smoothing operation on the target edge.

[0100] Strong constraints on radiotherapy-specific density levels refer to using a hard threshold corresponding to the density level as the constraint coefficient. In this embodiment, the smooth constraint coefficient for the density level corresponding to the air / cavity layer is: The smoothing constraint coefficient of the density level corresponding to the low-contrast soft tissue layer is The smoothing constraint coefficient of the density level corresponding to the tumor target layer is: The smoothing constraint coefficient of the density level corresponding to the bone layer is The smoothing constraint coefficient of the density level corresponding to the metal layer is .

[0101] In this embodiment, an artifact detection gradient threshold G is set. th =60 HU / voxel, calculate the three-dimensional gradient value and combine it with the density layer matrix to highlight the artifact residue areas around the simulated tumor target area and the brainstem and spinal cord, and generate the artifact mask matrix ArtifactMask; use 3×3×3 three-dimensional neighborhood mean smoothing, and perform secondary correction according to the radiotherapy-specific smoothing constraint coefficient γ (0.3 for the simulated tumor target area layer and 0.7 for the soft tissue layer where the brainstem / spinal cord is located) to obtain optimized correction data, and completely remove the artifacts around the target area and organs at risk.

[0102] Step S5: Output optimized correction data.

[0103] Step S5 specifically includes: converting the optimized correction data to a radiotherapy-specific format, outputting it in a clinically compatible format for tumor radiotherapy (such as DICOM-RT format). The output optimized correction data can be directly used for subsequent operations such as target delineation, positioning verification, image-guided localization, and dose calculation in tumor radiotherapy, completing the entire beam hardening correction process. Furthermore, the image of the optimized correction data as the output result can be directly integrated into the image guidance system of the radiotherapy linear accelerator without secondary format conversion.

[0104] In this embodiment, FinalData is converted to DICOM-RT format and directly output to the image guidance and target delineation system of the radiotherapy linear accelerator. After correction, the image is free of artifacts, the CT value is accurate, and the boundary between the simulated tumor target area and organs at risk is clear. It can be directly used for positioning verification, accurate target delineation and dose calculation in radiotherapy for head tumors.

[0105] The core parameters of this invention are designed according to the specific needs of image-guided localization in tumor treatment. They can be flexibly adjusted for different tumor sites such as lung cancer, head and neck tumors, and pelvic tumors, taking into account both radiotherapy accuracy and clinical versatility. The default values ​​and functions of the core parameters are shown in Table 1.

[0106] Table 1: Default values ​​and functions of core parameters

[0107]

[0108] This invention effectively addresses the core challenges of existing CBCT beam hardening correction methods in image-guided localization for tumor treatment through a three-tiered combined design: radiotherapy-specific three-dimensional density stratification, adaptive path correction, and precise optimization of target / organ-at-risk artifacts. It offers the following significant advantages:

[0109] 1. Improved Image Uniformity: Using a Catphan 700 phantom to test image quality uniformity (registration requirement standard is below 50 HU), the uniformity was substandard before applying the method of this invention, but far exceeded the registration standard after application. Comparative images are shown below. Figure 2A and Figure 2B As shown. By Figure 2A and Figure 2B It can be seen that the uniformity before correction is 67.997 HU, and the uniformity after correction is 7.968 HU (Note: the uniformity calculation method is Mean). max - Mean MIN ).

[0110] 2. Excellent artifact removal and strong 3D adaptability: Fully utilizing the interlayer correlation of 3D volumetric data, the residual artifact rate in the target area is reduced by more than 80%, adapting to tumors in multiple locations. Taking the head phantom as an example, the comparison images under the same window width and window level are as follows: Figures 3A-3D As shown. By Figures 3A-3D It can be seen that the present invention has a very good artifact removal effect.

[0111] 3. Highly efficient and real-time operation: 512×512×512 data processing time ≤3s, suitable for real-time guidance at the radiotherapy bedside;

[0112] 4. Easy integration and compatibility: It can be directly embedded into the radiotherapy CBCT post-processing system, outputting DICOM-RT format and seamlessly connecting to the radiotherapy workflow;

[0113] 5. Flexible and easy-to-use parameters: Default parameters are adapted to clinical practice and can be adjusted without professional algorithm knowledge.

[0114] The above description is merely a preferred embodiment of the present invention and is not intended to limit the scope of the invention. Various variations can be made to the above embodiments of the present invention. All simple and equivalent changes and modifications made in accordance with the claims and description of this application fall within the protection scope of the claims of this patent. All aspects not described in detail in this invention are conventional technical content.

[0115] The three-dimensional image post-processing method for CBCT beam hardening correction of the present invention constructs a synergistic system of CBCT beam hardening correction and detail preservation adapted to the image-guided localization requirements of tumor treatment, breaking through the limitations of the prior art, specifically in four aspects:

[0116] I. Dual Innovation in Processing Dimensions and Clinical Scenarios: The entire process of beam hardening correction is based on CBCT 3D volumetric data for tumor radiotherapy rather than 2D data. At the same time, a dedicated density stratification, correction model and parameter system are designed for the localization needs of tumor radiotherapy. Starting from the requirements of precision and real-time operation in tumor radiotherapy, and combining the characteristics of CBCT 3D volumetric data, a dedicated density stratification, correction model and full-process collaborative system are designed. This is different from existing general correction methods, and the processing logic is highly consistent with the requirements of precision and real-time operation in tumor radiotherapy.

[0117] II. Innovative Mechanism: This invention employs adaptive correction based on radiotherapy-specific three-dimensional voxel density stratification and path cumulative density analysis. A new tumor target layer is added and protected during correction. Combining the clinical needs of image-guided localization in tumor treatment with the physical characteristics of CBCT beam hardening, these methods are customized and modified specifically for tumor radiotherapy. Exclusive density stratification thresholds, path weights, correction models, and artifact optimization rules are designed, avoiding the limitations of traditional holistic correction and fixed-parameter correction. Compared to existing methods, artifact removal is cleaner, CT values ​​are more accurate, tumor details are better preserved, and the process is computationally efficient, easily integrated, and adaptable to the clinical workflow of tumor radiotherapy.

[0118] III. Innovative Process Design: A three-level, three-dimensional cascade closed loop is constructed, consisting of "radiotherapy-specific correction - target area / organ at risk artifact detection - radiotherapy-strong constraint optimization". This organically combines beam hardening artifact removal, CT value restoration, and preservation of details of low-contrast tumor lesions / organs at risk. It is also compatible with the post-processing system of tumor radiotherapy, forming an integrated CBCT image post-processing solution specifically for tumor radiotherapy, which solves the problem of the disconnect between existing technology and radiotherapy process.

[0119] IV. Artifact Optimization Innovation: A secondary artifact optimization method based on three-dimensional gradient detection and strong constraints of radiotherapy density layers is proposed. It focuses on the accurate identification and lightweight correction of artifacts around the tumor target area and organs at risk. By minimizing the smoothing operation of the target area layer, a balance between artifact removal and target area detail preservation is achieved, which is different from the existing indiscriminate global smoothing artifact removal methods.

[0120] In summary, the adaptive beam hardening correction method based on three-dimensional voxel density stratification of the present invention adopts a three-level three-dimensional cascade process of "three-dimensional density stratification and beam path analysis → radiotherapy-specific adaptive stratification correction → secondary optimization of target area / organs at risk artifacts". It removes artifacts and restores CT values ​​while preserving tumor and organ at risk details. It can be directly integrated into radiotherapy CBCT systems to meet the needs of precise radiotherapy for tumors.

[0121] The above description is merely a preferred embodiment of the present invention and is not intended to limit the scope of the invention. Various variations can be made to the above embodiments of the present invention. All simple and equivalent changes and modifications made based on the claims and description of this application fall within the protection scope of the claims of this patent. All aspects not described in detail in this invention are conventional technical content.

Claims

1. A three-dimensional image post-processing method for CBCT beam hardening correction, characterized in that, include: The CBCT 3D reconstructed volume data is subjected to 3D voxel density layering to generate a 3D density layering matrix, and the beam path cumulative density of each voxel is calculated. Based on the three-dimensional density layer matrix and path cumulative density, adaptive correction coefficients are calculated for CBCT three-dimensional reconstruction data, and adaptive correction of CT values ​​is performed to obtain preliminary correction data. The preliminary correction data were subjected to artifact residue detection and secondary optimization in the tumor target layer and low-contrast soft tissue layer to obtain optimized correction data.

2. The method according to claim 1, characterized in that, The three-dimensional density layer matrix DensityLayer includes a density layer label for each voxel, which is determined based on the CT value physical characteristics of key tissues in tumor radiotherapy and the clinical needs of tumor radiotherapy localization.

3. The method according to claim 2, characterized in that, Based on the physical characteristics of CT values ​​of key tissues in tumor radiotherapy, three-dimensional voxels are divided into five density levels, including air / cavity layer, low-contrast soft tissue layer, tumor target area layer, bone layer, and metal layer. The tumor target area layer is used to meet the precise positioning requirements of tumor radiotherapy. The CT values ​​for each density level are as follows: Air / cavity layer: CT < 500 HU; Low-contrast soft tissue layer: -500 HU < CT ≤ 200 HU; Tumor target area layer: 200 HU < CT ≤ 300 HU; Bone layer: 300 HU < CT ≤ 2000 HU; Metallic layer: 2000 HU < CT.

4. The method according to claim 1, characterized in that, The path cumulative density for each voxel is calculated, specifically including: based on the X-ray source location, detector parameters, and patient positioning from the CBCT scan, calculating the X-ray incident path for each voxel P(xx,yy,zz), and then calculating the path cumulative density D. sum (xx,yy,zz); Path cumulative density D sum The formula for calculating (xx,yy,zz) is: ， in, Let the coordinates be those of a voxel. The total number of voxels along the X-ray incident path, CT i The first X-ray in the incident path Raw CT values ​​of individual pixels, w i For the first The density hierarchy weighting coefficient of a voxel, where i is the ordinal number of the voxel.

5. The method according to claim 1, characterized in that, The larger the adaptive correction coefficient, the smaller the correction degree; the adaptive correction coefficient is related to the density level, and the adaptive correction coefficient corresponding to the density level of the tumor target area layer and the low contrast soft tissue layer is higher than the adaptive correction coefficient of other density levels; in addition, the adaptive correction coefficient is negatively correlated with the path cumulative density.

6. The method according to claim 5, characterized in that, The adaptive correction coefficients k(xx,yy,zz) are: ， Where k0(L) is the baseline correction coefficient for the density level, and L is the density level label. D is the attenuation adjustment coefficient. max This represents the maximum value of the path cumulative density. Let the coordinates be those of a voxel. For path cumulative density; Preliminary correction data for: ， in, These are the original CT values. For preliminary data correction, β is the adaptive correction coefficient, and β is the CT value offset compensation amount.

7. The method according to claim 6, characterized in that, When the density level label is air / cavity layer, k0(L) is When the density level label is a low-contrast soft tissue layer, k0(L) is... When the density level label is the tumor target layer, k0(L) is... When the density level label is the bone layer, k0(L) is... When the density level label is a metal layer, k0(L) is... .

8. The method according to claim 7, characterized in that, The preliminary corrected data underwent artifact residue detection and secondary optimization in the tumor target layer and low-contrast soft tissue layer to obtain optimized corrected data, specifically including: Based on the preliminary correction data, calculate the three-dimensional gradient value and detect and identify artifact regions based on the three-dimensional gradient value; For the artifact residue regions marked in the artifact mask matrix, fine correction is performed using three-dimensional neighborhood mean smoothing and strong constraints of radiotherapy-specific density levels to obtain optimized correction data.

9. The method according to claim 8, characterized in that, Based on the preliminary correction data, the three-dimensional gradient value is calculated, and artifact regions are detected and identified based on the three-dimensional gradient value. Specifically, this includes: A 3D neighborhood window is set for each voxel P(xx,yy,zz), and an artifact detection gradient threshold G is pre-set. th ; Based on the preliminary correction data and the three-dimensional neighborhood window, the three-dimensional gradient value G(xx,yy,zz) of the voxel is calculated to characterize the gray-level change rate between the voxel and its neighborhood. The three-dimensional gradient value G(xx,yy,zz) is compared with the artifact detection gradient threshold G. th The comparison is performed, and the artifact residue area is determined based on the comparison results; the artifact residue areas around the tumor target layer and the low-contrast soft tissue layer in the three-dimensional density layering matrix are marked to form an artifact mask matrix; Optimize calibration data for: ， Where N is the 3×3×3 effective neighborhood of the voxel being corrected, and M is the number of voxels in the effective neighborhood. It is a three-dimensional artifact mask matrix. For the initial data correction, γ(L) is the smoothing constraint coefficient of the density level, Q is the voxel, L is the density level label, and (xx,yy,zz) are the voxel coordinates.

10. The method according to claim 8, characterized in that, After obtaining the optimized correction data, the process also includes outputting the optimized correction data, specifically: converting the optimized correction data into a radiotherapy-specific format and outputting it in a format suitable for clinical use in tumor radiotherapy.

Images