Method and system for three-dimensional reconstruction of CT images based on attention mechanism
Patent Information
- Application Number
- CN202611064280.1
- Authority / Receiving Office
- CN · China
- Patent Type
- Patents(China)
- Current Assignee / Owner
- Filing Date
- 2026-07-17
- Publication Date
- 2026-10-09
- Estimated Expiration
- 2046-07-17
AI Technical Summary
[0002]在基于低剂量螺旋CT(LDCT)数据的肺结节三维重建任务中,现有方法大多采用逐层分割后堆叠或传统插值算法(如三线性插值、三次样条插值)生成三维体素网格;然而,由于低剂量采集协议下原始医学数字影像数据中往往存在较为明显的非均匀泊松噪声与硬化束伪影,同一结节区域在相邻切片上的灰度分布可能出现随机波动;在此类信噪比较低的场景下,传统体素重建手段有时难以稳定地捕捉结节边缘沿Z轴方向的微弱连续性特征;具体而言,当结节边缘的真实灰度梯度幅度与背景噪声梯度幅度较为接近时,现有算法容易将部分噪声波动误判为边缘结构,或者将原本连续的弱边缘信号平滑为渐变背景,从而在重建的三维模型表面引发非生理性的局部突变、细小孔洞或拓扑连接错误
[0022] The method efficiently parses the initial projection feature tensor using a multi-head self-attention encoding network, quickly capturing key features (edges, details) in the projection data and avoiding redundant computation of invalid features. The construction of the 3D volumetric attention mask, based on adaptive adjustment of local neighborhood features, makes the smoothing and compensation process more targeted, reducing unnecessary iterations and improving reconstruction accuracy while maintaining efficiency. In addition, the method introduces steps such as global numerical rate of change convergence judgment and dynamic range boundary truncation, which further enhance the stability of the reconstruction process and avoid problems such as iterative divergence and result fluctuations. The 2D spatial attention weight map generated by the multi-head self-attention encoding network can adaptively focus on the gradient features of the nodule edge region, strengthening the correlation of edge signals between adjacent slices. The constructed 3D volumetric attention mask can combine 3D local neighborhood features to perform anisotropic diffusion smoothing and detail compensation on the initially reconstructed voxel matrix, avoiding the smoothing of weak edge signals into a gradient background and preventing local abrupt changes, small holes, or topological connection errors caused by noise.
Smart Images

Figure CN122574164B_ABST
Abstract
Description
Technical Field
[0001] This invention relates to the field of data processing technology, and in particular to a method and system for three-dimensional reconstruction of CT images based on attention mechanisms. Background Technology
[0002] In the task of 3D reconstruction of lung nodules based on low-dose spiral CT (LDCT) data, most existing methods use layer-by-layer segmentation followed by stacking or traditional interpolation algorithms (such as trilinear interpolation and cubic spline interpolation) to generate a 3D voxel mesh. However, due to the presence of significant non-uniform Poisson noise and hardening bundle artifacts in the original medical digital image data under low-dose acquisition protocols, the gray-level distribution of the same nodule region on adjacent slices may exhibit random fluctuations. In such scenarios with low signal-to-noise ratios, traditional voxel reconstruction methods sometimes struggle to reliably capture the weak continuity features of the nodule edge along the Z-axis. Specifically, when the true gray-level gradient amplitude of the nodule edge is close to the background noise gradient amplitude, existing algorithms are prone to misjudging some noise fluctuations as edge structures or smoothing the originally continuous weak edge signals into a gradient background, thereby inducing non-physiological local mutations, small holes, or topological connection errors on the surface of the reconstructed 3D model.
[0003] For example, for a ground-glass nodule with a diameter of about 4 mm, the gray-scale gradient at its edge is relatively gentle. Under low-dose scanning conditions, the gray-scale value may fluctuate randomly by 10 to 20 Henle units between adjacent slices due to noise interference. When using conventional trilinear interpolation for voxel filling in the Z-axis direction, the interpolation result may mistakenly treat the noise peak as the nodule boundary and extend it outward by one pixel in a certain layer, while the real edge of the area in the next layer is treated as noise and suppressed. The surface of the reconstructed 3D model will show a non-realistic structure similar to a sawtooth or step. This geometric distortion may cause the calculated nodule volume to deviate from the true value in clinical monitoring, and also affect the reliable restoration of the spatial contact relationship between the nodule and the adjacent microvessels. Summary of the Invention
[0004] The technical problem to be solved by the present invention is to provide a method and system for three-dimensional reconstruction of CT images based on attention mechanism, which can improve the ability to preserve edge details and the overall image quality.
[0005] To solve the above-mentioned technical problems, the technical solution of the present invention is as follows:
[0006] Firstly, a CT image three-dimensional reconstruction method based on an attention mechanism, the method comprising:
[0007] Step 1: Collect sparse view sinusoidal data, analyze and standardize the sparse view sinusoidal data, and construct an initial projection feature tensor containing high-frequency edge gradient information and low-frequency background baseline information.
[0008] Step 2: Input the initial projection feature tensor into the pre-constructed multi-head self-attention encoding network to obtain the attention query feature tensor and the attention key feature tensor; calculate the dot product similarity between the transposes of the attention query feature tensor and the attention key feature tensor to obtain the two-dimensional attention weight map.
[0009] Step 3: Based on the two-dimensional spatial attention weight map, the initial projection feature tensor is subjected to element-wise weighted modulation processing to obtain the enhanced projection data stream;
[0010] Step 4: Substitute the enhanced projection data stream into the algebraic iterative reconstruction model to obtain the preliminary three-dimensional voxel matrix;
[0011] Step 5: Extract the three-dimensional local neighborhood features of the preliminary three-dimensional voxel matrix and construct a three-dimensional volumetric attention mask; use the three-dimensional volumetric attention mask to perform anisotropic diffusion smoothing and detail compensation on the preliminary three-dimensional voxel matrix to obtain the three-dimensional reconstructed image.
[0012] Secondly, attention-based CT image 3D reconstruction systems include:
[0013] The parsing module is used to collect sparse view sinusoidal data, parse and normalize the sparse view sinusoidal data, and construct an initial projection feature tensor containing high-frequency edge gradient information and low-frequency background baseline information.
[0014] The computation module is used to input the initial projected feature tensor into a pre-constructed multi-head self-attention encoding network to obtain the attention query feature tensor and the attention key feature tensor; and to calculate the dot product similarity between the transposes of the attention query feature tensor and the attention key feature tensor to obtain a two-dimensional spatial attention weight map.
[0015] The processing module is used to perform element-wise weighted modulation processing on the initial projection feature tensor based on the two-dimensional spatial attention weight map to obtain the enhanced projection data stream;
[0016] The reconstruction module is used to substitute the enhanced projection data stream into the algebraic iterative reconstruction model to obtain a preliminary three-dimensional voxel matrix.
[0017] The output module is used to extract the three-dimensional local neighborhood features of the preliminary three-dimensional voxel matrix and construct a three-dimensional volumetric attention mask. The three-dimensional volumetric attention mask is used to perform anisotropic diffusion smoothing and detail compensation on the preliminary three-dimensional voxel matrix to obtain a three-dimensional reconstructed image.
[0018] Thirdly, a computing device includes:
[0019] One or more processors;
[0020] A storage device for storing one or more programs that, when executed by one or more processors, cause the one or more processors to implement the method.
[0021] The above-described solution of the present invention has at least the following beneficial effects:
[0022] The method efficiently parses the initial projection feature tensor using a multi-head self-attention encoding network, quickly capturing key features (edges, details) in the projection data and avoiding redundant computation of invalid features. The construction of the 3D volumetric attention mask, based on adaptive adjustment of local neighborhood features, makes the smoothing and compensation process more targeted, reducing unnecessary iterations and improving reconstruction accuracy while maintaining efficiency. In addition, the method introduces steps such as global numerical rate of change convergence judgment and dynamic range boundary truncation, which further enhance the stability of the reconstruction process and avoid problems such as iterative divergence and result fluctuations. The 2D spatial attention weight map generated by the multi-head self-attention encoding network can adaptively focus on the gradient features of the nodule edge region, strengthening the correlation of edge signals between adjacent slices. The constructed 3D volumetric attention mask can combine 3D local neighborhood features to perform anisotropic diffusion smoothing and detail compensation on the initially reconstructed voxel matrix, avoiding the smoothing of weak edge signals into a gradient background and preventing local abrupt changes, small holes, or topological connection errors caused by noise. Attached Figure Description
[0023] Figure 1 This is a schematic diagram of the CT image three-dimensional reconstruction method based on the attention mechanism provided in the embodiments of the present invention.
[0024] Figure 2 This is a schematic diagram of an attention-based CT image three-dimensional reconstruction system provided in an embodiment of the present invention. Detailed Implementation
[0025] Exemplary embodiments of the present disclosure will now be described in more detail with reference to the accompanying drawings. While exemplary embodiments of the present disclosure are shown in the drawings, it should be understood that the present disclosure may be implemented in various forms and should not be limited to the embodiments set forth herein. Rather, embodiments are provided so that this disclosure will be thorough and complete, and will fully convey the scope of the disclosure to those skilled in the art.
[0026] like Figure 1 As shown, embodiments of the present invention propose a 3D reconstruction method for CT images based on an attention mechanism, the method comprising the following steps:
[0027] Step 1: Collect sparse view sinusoidal data, analyze and standardize the sparse view sinusoidal data, and construct an initial projection feature tensor containing high-frequency edge gradient information and low-frequency background baseline information.
[0028] Step 2: Input the initial projection feature tensor into the pre-constructed multi-head self-attention encoding network to obtain the attention query feature tensor and the attention key feature tensor; calculate the dot product similarity between the transposes of the attention query feature tensor and the attention key feature tensor to obtain the two-dimensional attention weight map.
[0029] Step 3: Based on the two-dimensional spatial attention weight map, the initial projection feature tensor is subjected to element-wise weighted modulation processing to obtain the enhanced projection data stream;
[0030] Step 4: Substitute the enhanced projection data stream into the algebraic iterative reconstruction model to obtain the preliminary three-dimensional voxel matrix;
[0031] Step 5: Extract the three-dimensional local neighborhood features of the preliminary three-dimensional voxel matrix and construct a three-dimensional volumetric attention mask; use the three-dimensional volumetric attention mask to perform anisotropic diffusion smoothing and detail compensation on the preliminary three-dimensional voxel matrix to obtain the three-dimensional reconstructed image.
[0032] In this embodiment of the invention, the initial projection feature tensor is efficiently parsed using a multi-head self-attention coding network, which can quickly capture key features (edges, details) in the projection data and avoid redundant calculations of invalid features. The construction of the 3D volumetric attention mask is based on adaptive adjustment of local neighborhood features, making the smoothing and compensation process more targeted, reducing unnecessary iterations, and improving reconstruction accuracy while taking into account reconstruction efficiency. In addition, the steps introduced in the method, such as global numerical rate of change convergence judgment and dynamic range boundary interception, further improve the stability of the reconstruction process and avoid problems such as iterative divergence and result fluctuations. The 2D spatial attention weight map generated by the multi-head self-attention coding network can adaptively focus on the gradient features of the nodule edge region and strengthen the correlation of edge signals between adjacent slices. The constructed 3D volumetric attention mask can combine 3D local neighborhood features to perform anisotropic diffusion smoothing and detail compensation on the initially reconstructed voxel matrix, which not only avoids weak edge signals being smoothed into a gradient background, but also prevents local abrupt changes, small holes, or topological connection errors caused by noise.
[0033] In a preferred embodiment of the present invention, step 1, acquiring sparse viewpoint sinusoidal image data, analyzing and standardizing the sparse viewpoint sinusoidal image data, and constructing an initial projection feature tensor containing high-frequency edge gradient information and low-frequency background baseline information, may include:
[0034] Step 101: Acquire sparse-view sinusoidal data. Perform logarithmic domain transformation and channel consistency compensation on the sparse-view sinusoidal data to obtain linearized corrected projection data. Specifically, this includes: starting a low-dose spiral CT scanner and performing sparse-view scanning acquisition on the lung region. During the scanning process, the scan dose is strictly controlled within the clinical low-dose standard range, typically controlled at 30% to 50% of the conventional CT scan dose. This ensures that the acquired sparse-view sinusoidal data can completely preserve the effective imaging information of normal lung tissue and small lesions such as ground-glass nodules, while minimizing the radiation dose harm to the human body. The acquired sparse-view sinusoidal data is stored in the form of a two-dimensional matrix, with the matrix rows... For different detection channels, columns correspond to different detection angles. Each matrix element represents the signal value of X-rays penetrating lung tissue at the corresponding detection channel and angle. Due to the low-dose acquisition characteristics, this data inevitably carries non-uniform Poisson noise and hardening beam artifacts. These interferences cause irregular fluctuations in signal values, affecting the accuracy of subsequent feature extraction and hindering the precise capture of edge features of small lesions. After acquisition, a logarithmic domain transformation is performed on the sparse viewpoint sine curve data. Specifically, for each sampling point in the sine curve, its original signal value is extracted one by one, and the natural logarithm of each original signal value is taken individually. The specific calculation process is as follows: Let the original signal value of a certain sampling point be... The value obtained by taking the natural logarithm is The calculation process has a defined expression, which is: =ln( ),in This represents the raw X-ray penetration signal value at each sampling point in the sparse view sinogram, without any processing. The linearized signal value obtained after the original signal value is transformed by natural logarithm is represented. This processing method can effectively eliminate the exponential characteristics of X-ray attenuation in lung tissue, convert the originally nonlinear attenuation signal into a linear signal, making the signal change pattern easier to capture, and facilitating the subsequent extraction of high-frequency edge gradient information and low-frequency background baseline information.
[0035] After the logarithmic domain transformation, preliminary data analysis revealed significant signal strength differences between different detection channels. These differences stem from variations in the hardware performance of the detection channels themselves, leading to inconsistent signal values for the same area collected by different channels. This affects the consistency of subsequent feature extraction, necessitating channel consistency compensation. The specific procedure is as follows: First, statistically analyze the signal data from all detection channels. Assuming there are several detection channels, each with several sampling points, calculate the average value of all sampling points within each channel to obtain the individual signal mean for each channel. Specifically, for any detection channel, its individual signal mean equals the sum of all sampling point values within that channel divided by the number of sampling points in that channel. That is, the individual signal mean equals the sum of the first and second sampling point values of that channel, and so on, accumulating up to the last sampling point value. Then, divide the sum by the total number of sampling points in that channel. Repeat this process to calculate the individual signal mean for all detection channels. The process involves: 1) Calculating the sum of the signal mean values of all detection channels, i.e., summing the individual signal mean values of all detection channels sequentially; 2) Dividing this sum by the total number of detection channels to obtain the overall signal mean value of all detection channels, i.e., the overall signal mean value equals the sum of the mean values divided by the total number of detection channels; 3) For each detection channel, calculating the difference between the individual signal mean value of that channel and the overall signal mean value of all channels to obtain the compensation difference value for each channel. Specifically, the compensation difference value for any detection channel equals the overall signal mean value minus the individual signal mean value of that channel; 4) Adding the corresponding compensation difference value to the values of all sampling points within that channel. Specifically, the compensated value of any sampling point in that channel equals the original value of that sampling point plus the compensation difference value for that channel. This process is strictly performed for each detection channel to ensure that the signal strength of all detection channels remains consistent, completely eliminating errors between channels; 5) After the above two steps of logarithmic domain transformation and channel consistency compensation, the linearized corrected projection data is finally obtained.
[0036] Step 102: Based on the linearized corrected projection data, extract the neighborhood numerical difference of each sampling point along the preset detection direction. Perform gradient amplitude aggregation processing on the neighborhood numerical difference to obtain high-frequency edge gradient information. Specifically, this includes: determining the preset detection direction. Combining the detector array layout of low-dose spiral CT scanning, the preset detection direction is set to be completely consistent with the detector array direction, specifically covering the horizontal and vertical detection directions. The two directions are perpendicular to each other. The horizontal detection direction is consistent with the arrangement direction of the detector channels, and the vertical detection direction is consistent with the direction of change of the detection angle, ensuring that the surrounding area of each sampling point can be fully covered and the edges in different directions in the image can be accurately captured. Edge features are used to avoid missing edge features. After determining the preset detection direction, for each sampling point in the linearized corrected projection data, neighboring sampling points are selected in both the horizontal and vertical detection directions. The number of neighboring sampling points is flexibly determined according to the sampling density of the sparse view. When the sampling density is low, the number of neighboring sampling points is appropriately increased, usually two on the left and two on the right, for a total of four neighboring sampling points. When the sampling density is high, the number of neighboring sampling points is appropriately reduced, usually one on the left and one on the right, for a total of two neighboring sampling points, to ensure that the selected neighboring sampling points can accurately reflect the signal changes around the sampling point. After selection, the numerical difference between each sampling point and each neighboring sampling point is calculated. The calculation method is uniformly the linearized corrected projection data value of the sampling point minus the corresponding neighboring sampling point value. The specific calculation process is as follows: let the value of the current sampling point be the current sampling value, and the value of a certain neighboring sampling point be the neighboring sampling value. The numerical difference between the two is equal to the current sampling value minus the neighboring sampling value. Each sampling point will obtain multiple neighboring numerical differences in the horizontal detection direction and multiple neighboring numerical differences in the vertical detection direction. All these differences are collected as the initial gradient correlation data of the sampling point.
[0037] Gradient magnitude aggregation is performed on all neighborhood numerical differences for each sampling point. The specific process involves calculating the absolute value of each neighborhood numerical difference individually to eliminate the positive and negative influences of the difference, retaining only its magnitude. Specifically, a neighborhood numerical difference is defined as the difference data, and its absolute value is the absolute value of the difference data. The absolute values of all neighborhood numerical differences for that sampling point are summed, and the sum is the gradient magnitude for that sampling point. Specifically, if a sampling point has several neighborhood numerical differences, its gradient magnitude is equal to the absolute value of the first neighborhood numerical difference plus the absolute value of the second neighborhood numerical difference, and so on until the absolute value of the last neighborhood numerical difference. The magnitude of the gradient magnitude reflects the intensity of signal changes around the sampling point; the more intense the signal change, the larger the gradient magnitude, and the more obvious the edge features in the corresponding image. Conversely, the more gradual the signal change, the smaller the gradient magnitude, corresponding to background areas in the image. The gradient magnitudes of all sampling points in the linearized corrected projection data are then combined according to their original spatial locations to form high-frequency edge gradient information.
[0038] Step 103 involves subtracting the corresponding high-frequency edge gradient information values from the sampling point values in the linearized corrected projection data to obtain the low-frequency background baseline information. Specifically, this includes: clarifying the spatial correspondence between the linearized corrected projection data and the high-frequency edge gradient information. Both are two-dimensional matrices with the same size as the original sparse viewpoint sine wave data, with rows and columns corresponding to the same detection channels and detection angles, ensuring that the spatial coordinates of each sampling point in the two datasets are completely consistent, avoiding calculation errors caused by positional misalignment; for each sampling point in the linearized corrected projection data, finding the gradient in the high-frequency edge gradient information that completely corresponds to the spatial position of that sampling point. The gradient magnitude value, i.e., the gradient magnitude value of a sampling point in a certain row and column of the linearized corrected projection data, corresponds to the gradient magnitude value in the same row and column of the high-frequency edge gradient information, ensuring that each sampling point can be matched with a unique gradient magnitude value. A subtraction operation is performed. The specific calculation process is as follows: let the value of a sampling point in the linearized corrected projection data be the corrected sampled value, and the gradient magnitude value at its corresponding position be the gradient magnitude. The value obtained after subtraction is the background sampled value, that is, the background sampled value equals the corrected sampled value minus the gradient magnitude. The above subtraction operation is performed individually on each sampling point, without omitting any sampling point. After subtraction, the subtracted value of each sampling point is obtained. The subtracted values of all sampling points are combined according to their original spatial positions to form low-frequency background baseline information. The core function of this process is to effectively separate the high-frequency edge features and low-frequency background features in the linearized corrected projection data. The low-frequency background baseline information mainly reflects the overall background signal of lung tissue, including the basic signal of normal lung tissue.
[0039] Step 104: Based on the high-frequency edge gradient information, extract the global maximum and global minimum values. Map the values of each sampling point in the high-frequency edge gradient information to a preset standard interval according to the range ratio to obtain standardized high-frequency edge gradient information. Specifically, this includes: performing a full traversal of all sampling point values in the high-frequency edge gradient information, comparing the gradient magnitude value of each sampling point one by one, selecting the sampling point value with the largest value and determining this value as the global maximum value of the high-frequency edge gradient information, and simultaneously selecting the sampling point value with the smallest value and determining this value as the global minimum value of the high-frequency edge gradient information; after the traversal is completed, calculate the global maximum and global minimum values. The difference between the global minimum and the global maximum is the range of the high-frequency edge gradient information. Specifically, the range is calculated as the global maximum minus the global minimum. The range reflects the overall fluctuation range of the values in the high-frequency edge gradient information. A preset standard interval is set according to the conventional needs of CT image processing and the input requirements of the subsequent multi-head self-attention coding network. It is usually set to [0,1] to ensure that the standardized values can fully adapt to the network's input range, avoiding feature loss or processing distortion due to values exceeding the network's processing range. For each sampling point value in the high-frequency edge gradient information, standardization calculations are performed separately. The specific calculation process is as follows: The difference between the sampled point value and the global minimum is calculated to obtain the offset of the sampled point relative to the global minimum. Specifically, the gradient magnitude of a sampled point is taken as the gradient value, and its offset is equal to the gradient value minus the global minimum. This offset is then divided by the range of the high-frequency edge gradient information to obtain the proportionality coefficient corresponding to the sampled point value. The proportionality coefficient is calculated as the offset divided by the range. The proportionality coefficient reflects the relative position of the sampled point value within the overall fluctuation range, and its value ranges from [0,1]. Finally, this proportionality coefficient is multiplied by the length of a preset standard interval to obtain the relative position of the sampled point value within the preset standard interval. The numerical values are calculated as follows: the length of the preset standard interval [0,1] is 1 - 0 = 1, and the relative value is equal to the scaling factor multiplied by 1. The minimum value of the preset standard interval is added to the relative value to obtain the standardized value of the sampling point. The standardized value is equal to the relative value plus 0, which is equal to the scaling factor. That is, the standardized value is equal to the difference between the gradient value and the global minimum value, and then divided by the range between the global maximum value and the global minimum value. After all sampling points have completed the above standardization calculation, the standardized values of all sampling points are combined according to their original spatial positions to form standardized high-frequency edge gradient information.
[0040] Step 105: Based on the low-frequency background baseline information, extract the global maximum and global minimum values. Map the values of each sampling point in the low-frequency background baseline information to a preset standard interval according to the range ratio to obtain standardized low-frequency background baseline information. Concatenate and stitch the standardized high-frequency edge gradient information with the standardized low-frequency background baseline information to obtain the initial projection feature tensor. Specifically, this includes processing the low-frequency background baseline information according to the same specific operation method as in step 104. The specific operation process is as follows: comprehensively traverse all sampling point values in the low-frequency background baseline information, compare the values of each sampling point one by one, and select the global maximum and global minimum values; calculate the range between the global maximum and global minimum values. The specific calculation process is as follows: the range equals the global maximum value minus the global minimum value; for each sampling point value, the difference between this value and the global minimum value is calculated, specifically, this difference equals the low-frequency background baseline information value of this sampling point minus the global minimum value, where this value is the low-frequency background baseline information value of this sampling point; the difference is divided by the range to obtain the scaling factor, specifically, the scaling factor equals the difference divided by the range; the scaling factor is multiplied by the interval length of the preset standard interval, and added to the minimum value of the preset standard interval to obtain the standardized value of this sampling point, specifically, the preset standard interval is still [0,1], the interval length is 1, therefore the standardized value equals the scaling factor multiplied by 1 plus... 0, meaning the standardized value is equal to the difference between the low-frequency background baseline information value of the sampling point and the global minimum value, then divided by the range of the global maximum value minus the global minimum value; after all sampling points have completed the standardization calculation, standardized low-frequency background baseline information is formed, and the preset standard interval used in this standardization process is completely consistent with the preset standard interval of the high-frequency edge gradient information in step 104, both being [0,1], ensuring that the numerical intervals of high-frequency edge features and low-frequency background features remain consistent, avoiding numerical imbalance problems in the subsequent cascade stitching process, and ensuring that the stitched feature tensor can be stably input into the subsequent network for processing; after the standardization process is completed, the standardized high-frequency edge gradient information and the standardized low-frequency background... The baseline information is cascaded and stitched together. Specifically, the standardized low-frequency background baseline information is used as the first channel of the initial projection feature tensor, and the standardized high-frequency edge gradient information is used as the second channel of the initial projection feature tensor. During the stitching process, it is strictly ensured that the spatial dimensions of the two channels are completely consistent, that is, the number of rows and columns of the two channels are exactly the same, and the spatial positions of the sampling points in each channel correspond one-to-one without any misalignment. After stitching, an initial projection feature tensor containing high-frequency edge gradient information and low-frequency background baseline information is formed. The dimension of this tensor is the number of rows multiplied by the number of columns multiplied by two. The first channel stores the standardized low-frequency background baseline information, and the second channel stores the standardized high-frequency edge gradient information.
[0041] This embodiment separates and extracts high-frequency edge gradient information from low-frequency background baseline information, and uses numerical subtraction calculation to completely separate edge features from background signals. This can accurately distinguish lesion edge features from lung tissue background signals, effectively prevent edge features from being masked by background noise, and also prevent background signals from interfering with edge feature extraction.
[0042] In a preferred embodiment of the present invention, step 2, inputting the initial projection feature tensor into a pre-constructed multi-head self-attention encoding network, obtains an attention query feature tensor and an attention key feature tensor; calculating the dot product similarity between the transposes of the attention query feature tensor and the attention key feature tensor to obtain a two-dimensional spatial attention weight map, may include:
[0043] Step 201 involves inputting the initial projection feature tensor into a pre-constructed multi-head self-attention coding network. The initial projection feature tensor is then linearly mapped using the network's preset query weight parameters to obtain the initial query feature data. Specifically, this includes constructing the multi-head self-attention coding network, which is specifically designed for the initial projection feature tensor of low-dose CT images. This network is adapted to the characteristics of low-dose CT images, such as high noise levels and indistinct features of small lesions. The overall structure consists of an input layer, a linear mapping layer, a grouping and stitching layer, an axis interchange layer, a numerical computation layer, and an output layer connected sequentially. All layers are fully connected to ensure lossless feature transfer. The input dimension of the input layer strictly matches the dimension of the initial projection feature tensor. The initial projected feature tensor has a dimension of the number of rows multiplied by the number of columns multiplied by two. Therefore, the input dimension of the input layer is also set to the number of rows multiplied by the number of columns multiplied by two to ensure that the initial projected feature tensor can be completely input into the network without dimension misalignment or feature loss. The linear mapping layer contains two independent linear mapping branches, corresponding to the calls to query weight parameters and key weight parameters, respectively, and is used to generate initial query feature data and initial key feature data. The grouping and concatenation layer is used to group, split, and reassemble the feature data after linear mapping. The axis interchange layer is used to realize the dimension transformation of the feature tensor. The numerical operation layer is used to perform similarity calculation, scaling, exponential mapping, and normalization operations. The output layer is used to output the final two-dimensional spatial attention weight map.
[0044] After the network was built, training was immediately initiated. The training dataset specifically selected the initial projection feature tensor corresponding to the low-dose CT sparse-view sine wave map preprocessed in steps 101 to 105. It was also supplemented with manually annotated reference maps of attention weights at the edges of small lung lesions as annotation data. The annotation data was manually labeled according to the actual location, size, and edge features of the lung lesions to ensure accuracy. Gradient descent optimization was used during training, with the learning rate dynamically adjusted based on training iterations. The initial learning rate was set to 0.001, decreasing to 0.9 times the original rate every 100 iterations. The mean squared error between the predicted attention weights and the manually annotated data was used as the loss function. The loss function value was calculated after each iteration. If the loss function value was greater than a preset threshold, training continued; if the loss function value was less than or equal to the preset threshold, training stopped. At this point, network training was complete, and all network parameters were saved, including the preset query weight parameters and key weight parameters. The query weight parameters were three-dimensional matrices precisely matched to the dimensions of the initial projection feature tensor. The three dimensions of the matrix corresponded to the feature dimensions of the initial projection feature tensor, the multi-head auto-weights, and the multi-head auto-weights, respectively. The preset number of heads and feature dimensions of each head in the attention encoding network ensure that the linear mapping process can be executed smoothly without dimension mismatch issues. After the multi-head self-attention encoding network is deployed after training, the initial projection feature tensor obtained in step 105 is completely input into the input layer of the network. After input, the network automatically calls the preset query weight parameters saved during training to perform linear mapping processing on the initial projection feature tensor. Specifically, for each sampling point in the initial projection feature tensor, the feature value of the sampling point is extracted, and then the parameter value in the query weight parameters that completely corresponds to the spatial position and feature dimension of the sampling point is extracted. The feature value of the sampling point is multiplied by the corresponding parameter value to obtain a single multiplication result. All multiplication results corresponding to all feature dimensions of the sampling point are accumulated. The accumulated value is the initial query feature data corresponding to the sampling point. The above multiplication and accumulation operation is performed on each sampling point individually, without omitting any sampling point. The initial query feature data of all sampling points are combined together according to their original spatial positions to form complete initial query feature data.
[0045] Step 202: Perform grouping and segmentation on the initial query feature data, splitting it into multiple independent query feature segments; then, concatenate and reassemble these segments according to a preset data arrangement to obtain the attention query feature tensor. Specifically, this includes: determining the number of group segments, which is strictly consistent with the preset number of heads in the multi-head self-attention encoding network. The preset number of heads is set based on the complexity of the low-dose CT image features and the size and number of small lesions, typically set to 8 or 16 heads. A higher preset number of heads allows for more comprehensive feature dimensions to be captured, resulting in more accurate capture of the edge features of small lesions such as ground-glass nodules in the lungs. Considering the actual processing needs of low-dose CT images, a setting of 8 heads is preferred. During the grouping and segmentation process, the initial query feature data is uniformly segmented according to the feature dimensions. The total number of feature dimensions of the initial query feature data is calculated, and then divided by the preset number of heads to obtain the number of feature dimensions for each independent query feature segment, ensuring that the number of feature dimensions for each segment is exactly the same, achieving uniform segmentation. During the segmentation process, strict adherence to the preset number of heads is maintained for each independent query feature segment. The spatial dimension of the query feature fragments is kept consistent with the spatial dimension of the initial query feature data. That is, the number of rows and columns of each independent query feature fragment is exactly the same as the initial query feature data, only the number of feature dimensions differs, avoiding data loss, dimension disorder, or spatial position shift after segmentation. After segmentation, multiple independent query feature fragments with the same number of preset heads are obtained. Each independent query feature fragment corresponds to a unique group number, which starts from 1 and increases sequentially. The splicing and recombination operation is performed according to the preset data arrangement order, which is set to be in ascending order of group number. The specific operation of splicing and recombination is to sequentially superimpose the spatial dimensions of each independent query feature fragment. During the splicing process, the feature dimensions of each independent query feature fragment remain unchanged, and the feature values and spatial positions within each fragment are not changed. Only the spatial dimensions of each fragment are concatenated and spliced. After splicing, an attention query feature tensor is obtained. The dimension of this tensor is the preset number of heads multiplied by the number of rows in the spatial dimension of the initial query feature data multiplied by the number of columns in the spatial dimension of the initial query feature data multiplied by the feature dimension of each head.
[0046] Step 203: Input the initial projected feature tensor into the pre-constructed multi-head self-attention encoding network. Perform linear mapping on the initial projected feature tensor using the network's preset key weight parameters to obtain initial key feature data. Perform grouping and segmentation on the initial key feature data to split it into multiple independent key feature segments. Specifically, input the initial projected feature tensor back into the input layer of the trained and deployed multi-head self-attention encoding network. Unlike the operation in step 201, this time the network calls the preset key weight parameters saved during training, rather than querying the weight parameters. The key weight parameters are also three-dimensional matrices, with their three dimensions corresponding to the feature dimension of the initial projected feature tensor, the preset number of heads, and the feature dimension of each head, respectively. The dimensions are completely consistent with the dimensions of the query weight parameters, ensuring accurate adaptation with the dimensions of the initial projected feature tensor and avoiding the problem of dimension mismatch during linear mapping.
[0047] The linear mapping process is completely consistent with the linear mapping process for the initial query feature data in step 201. Specifically, for each sampling point in the initial projected feature tensor, the feature value of that sampling point is extracted, and then the parameter value in the key weight parameters that completely corresponds to the spatial position and feature dimension of that sampling point is extracted. The feature value of that sampling point is multiplied by the corresponding parameter value to obtain a single multiplication result. All multiplication results corresponding to all feature dimensions of that sampling point are accumulated. The accumulated value is the initial key feature data corresponding to that sampling point. The above multiplication and accumulation operation is performed on each sampling point individually. The initial key feature data of all sampling points are combined together according to their original spatial positions to form complete initial key feature data. The dimension of this initial key feature data is completely consistent with the dimension of the initial query feature data, that is, the number of rows and columns in the spatial dimension, as well as the total number of feature dimensions, are the same as those in the initial query feature data, which is used for subsequent... The grouping and reassembly operations provide convenience, ensuring that the grouped segments can accurately match the segments of the initial query feature data. After obtaining the initial key feature data, the grouping operation is performed in the same way as in step 202. The number of groups is also consistent with the preset number of heads in the multi-head self-attention encoding network. During the segmentation process, the initial key feature data is also uniformly segmented according to the feature dimensions. The total number of feature dimensions of the initial key feature data is calculated. The total number of feature dimensions is divided by the preset number of heads to obtain the number of feature dimensions of each independent key feature segment. This ensures that the number of feature dimensions of each independent key feature segment is the same and that the spatial dimension is consistent with the initial key feature data. After the segmentation is completed, multiple independent key feature segments with the same number of heads as the preset number are obtained. The dimension of each independent key feature segment is completely consistent with the dimension of the independent query feature segment obtained in step 202, that is, the number of rows and columns in the spatial dimension and the number of feature dimensions are all the same.
[0048] Step 204 involves concatenating and recombining multiple independent key feature fragments according to a preset data arrangement order to obtain an attention key feature tensor. Specifically, the preset data arrangement order is identical to that used in Step 202, both being ordered from smallest to largest group number. This ensures that the structure of the concatenated and recombined feature tensor remains consistent with that of the attention query feature tensor. The concatenation and recombining method is also identical to that used in Step 202. Specifically, the spatial dimensions of each independent key feature fragment are sequentially superimposed. During the concatenation process, the group number of each independent key feature fragment is checked one by one, concatenating sequentially according to number 1, number 2, and so on until the last number. After concatenation, the dimensions of the concatenated feature tensor are checked again to ensure that the dimensions of the concatenated attention key feature tensor are completely consistent with those of the attention query feature tensor, both being the preset number of heads multiplied by the number of rows in the spatial dimension of the initial projected feature tensor multiplied by the number of columns in the spatial dimension of the initial projected feature tensor multiplied by the feature dimension of each head. Through this concatenation and recombining operation, a structurally sound and dimensionally accurate attention key feature tensor is obtained.
[0049] Step 205 involves performing an axis swapping operation on the spatial dimension and feature dimension of the attention key feature tensor to obtain a transposed key feature tensor. Specifically, this includes clarifying the existing dimensional order of the attention key feature tensor. After the concatenation and recombination in step 204, the dimensional order of the attention key feature tensor is: preset number of heads, number of rows in the spatial dimension, number of columns in the spatial dimension, and feature dimension, arranged sequentially. The preset number of heads is the first dimension, the number of rows in the spatial dimension is the second dimension, the number of columns in the spatial dimension is the third dimension, and the feature dimension is the fourth dimension. Since the attention key feature tensor and the attention query feature tensor need to be used for dot product similarity calculation later, the existing dimensional order cannot be changed. To achieve accurate computation between the two, it is necessary to perform an axis swapping operation on the attention key feature tensor between the spatial dimension and the feature dimension. The specific operation of axis swapping is to swap the axis of the fourth dimension (the feature dimension) of the attention key feature tensor with the axis of the third dimension (the number of columns in the spatial dimension). After the swap, the feature values in each dimension are not changed, only the order of the four dimensions is changed. After the axis swap, the dimension order of the transposed key feature tensor becomes the preset number of heads, the number of rows in the spatial dimension, the feature dimension, and the number of columns in the spatial dimension. The preset number of heads is still the first dimension, the number of rows in the spatial dimension is still the second dimension, the feature dimension becomes the third dimension, and the number of columns in the spatial dimension becomes the fourth dimension.
[0050] Step 206 involves performing numerical multiplication and summation operations on the corresponding spatial locations of the attention query feature tensor and the transposed key feature tensor to obtain the original similarity data. Specifically, this includes: clarifying the dimensional correspondence between the attention query feature tensor and the transposed key feature tensor. The dimensional order of the attention query feature tensor is: preset number of heads, spatial dimension (rows), spatial dimension (columns), and feature dimension. The dimensional order of the transposed key feature tensor is: preset number of heads, spatial dimension (rows), feature dimension, and spatial dimension (columns). The preset number of heads and spatial dimension (rows) are completely identical, and the feature dimension and spatial dimension (columns) correspond to each other. This dimensional correspondence enables precise spatial location operations, ensuring accurate matching and computation of features at each spatial location. During the computation, numerical multiplication and summation operations are performed independently for each preset head, without mixed operations across preset heads, ensuring that the computation results for each preset head are independent and accurate. The specific computation process is as follows: first, for a specific... A preset head is defined, and a spatial location is extracted from the attention query feature tensor corresponding to that head. This spatial location is determined by the number of rows and columns in the spatial dimension. Then, all values in the transposed key feature tensor corresponding to that head that are identical to that spatial location and have the same feature dimension are extracted. Each feature value at that spatial location in the attention query feature tensor is multiplied one by one with each value of that spatial location and corresponding feature dimension in the transposed key feature tensor, resulting in multiple multiplication results, one for each feature dimension. The multiplication results corresponding to all feature dimensions at the same spatial location are summed, and the single sum is the original similarity value corresponding to that preset head and spatial location. The above multiplication and summation operation is performed separately for each preset head and each spatial location, without omitting any preset head or any spatial location. The original similarity values of all preset heads and all spatial locations are combined according to their corresponding positions to form the original similarity data.
[0051] Step 207: Perform a division operation with a preset scaling factor on the original similarity data to obtain scaled similarity data; then perform an exponential function mapping operation on the scaled similarity data to obtain positive exponential similarity data. Specifically, this includes: determining the preset scaling factor, which is strictly combined with the preset number of heads in the multi-head self-attention encoding network and the feature dimension of each head. The specific calculation method is the reciprocal of the square root of the feature dimension of each head, and has a clear expression: ,in This represents the preset proportional coefficient. This represents the feature dimension of each head. For example, if the feature dimension of each head is... =64, then the preset scaling factor of this head is 64. =0.125; The core purpose of setting this scaling factor is to avoid the original similarity data being too large, which could lead to numerical saturation in the subsequent exponential function mapping process, and thus avoid the problem of unclear feature similarity distinction. This aligns with the characteristics of low-dose CT images, where the feature value range is relatively small and precise feature similarity distinction is required. After the scaling factor is determined, each value in the original similarity data is divided independently without batch calculation to ensure calculation accuracy. Specifically, each original similarity value is divided individually by the preset scaling factor, which has a clear expression: ,in This represents the scaled similarity score. Represents the original similarity score. This represents a preset scaling factor. Each original similarity value corresponds to a scaled similarity value, avoiding many-to-one or one-to-many relationships. After division, all original similarity values are combined to form scaled similarity data. The dimensions of this scaled similarity data are completely consistent with those of the original similarity data; that is, the preset number of heads, spatial dimension rows, spatial dimension columns, and spatial dimension columns remain unchanged, only the values at each position are scaled. After scaling, an exponential function mapping process is performed on the scaled similarity data. The natural exponential function is used for mapping. Specifically, each value in the scaled similarity data is used individually as the exponent of the natural exponential function, and the result of the natural exponential function is calculated to obtain the positive exponential similarity value. The expression is: ,in This represents the positive exponential similarity score, where e represents the natural constant. This represents the scaled similarity value. For example, if a certain scaled similarity value... =0.5, then the corresponding positive exponential similarity value is Each scaled similarity value is processed separately, without missing any value. All positive exponential similarity values are combined together according to their original positions to form positive exponential similarity data.
[0052] Step 208: Accumulate the positive exponential similarity data position by position along the preset feature dimension direction to obtain the normalized denominator data corresponding to each spatial position; divide the value of each spatial position in the positive exponential similarity data by the corresponding normalized denominator data to obtain a two-dimensional spatial attention weight map. Specifically, this includes: determining the preset feature dimension direction. Combining the dimensional characteristics of the positive exponential similarity data, the preset feature dimension direction is set as the last spatial dimension direction in the positive exponential similarity data, that is, the direction of the number of columns in the spatial dimension. Perform position-by-position numerical accumulation along this direction. The specific process of position-by-position numerical accumulation is as follows: for each preset header, each number of rows in each spatial dimension, and each number of columns in each spatial dimension in the positive exponential similarity data, determine the spatial coordinates corresponding to the position, and then extract all positive exponential similarity values along the preset feature dimension direction, that is, the direction of the last number of columns in the spatial dimension, and accumulate these values one by one. The single value obtained after accumulation is the normalized denominator data corresponding to the preset header and the spatial position.
[0053] The above accumulation operation is performed individually for each preset head and each spatial location. All normalized denominator data are combined according to their corresponding positions to form a normalized denominator data set. The dimensions of this set are the number of preset heads, the number of rows in the spatial dimension, and the number of columns in the spatial dimension, consistent with the first three dimensions of the positive exponential similarity data. After obtaining the normalized denominator data, normalization processing is performed. Specifically, the value of each spatial location in the positive exponential similarity data is divided by the normalized denominator data corresponding to that spatial location. That is, for the same preset head, the same number of rows in the spatial dimension, and the same number of columns in the spatial dimension, the positive exponential similarity value is divided by the corresponding normalized denominator data to obtain... The result is the attention weight value corresponding to that spatial location. During the division process, the position of each positive exponential similarity value and the corresponding normalized denominator data are strictly checked to ensure a one-to-one correspondence and to avoid weight distortion caused by matching errors. At the same time, the case of a denominator of zero is avoided. If the special case of a denominator of zero occurs, the attention weight value of the corresponding position is set to 0 to ensure the integrity of the data. The attention weight values of all spatial locations are combined together according to their original spatial locations to form a two-dimensional spatial attention weight map. The dimension of this weight map is completely consistent with the spatial dimension of the initial projected feature tensor, that is, the number of rows and columns are the same as the initial projected feature tensor.
[0054] This embodiment obtains an accurate two-dimensional spatial attention weight map through a series of operations such as axis interchange, numerical multiplication and accumulation, scaling, exponential mapping and normalization. It can effectively distinguish between lesion edges, background areas and noise areas, highlight the importance of small lesion features and suppress the interference of noise and irrelevant background.
[0055] In a preferred embodiment of the present invention, step 3, which involves performing element-wise weighted modulation processing on the initial projection feature tensor based on a two-dimensional spatial attention weight map to obtain an enhanced projection data stream, may include:
[0056] Step 301 involves matching the two-dimensional spatial attention weight map according to the spatial position distribution of the initial projection feature tensor to obtain spatially aligned weight data. Specifically, this includes: clarifying the spatial dimension correspondence between the two-dimensional spatial attention weight map and the initial projection feature tensor, ensuring their spatial dimensions are identical (the number of rows and columns corresponding to the original sparse viewpoint sine wave data). The number of rows corresponds to the number of detection channels in low-dose spiral CT scanning, and the number of columns corresponds to the number of detection angles. Each weight value in the two-dimensional spatial attention weight map is extracted, and its corresponding spatial position is determined by the corresponding detection channel and angle. Then, a spatial position in the initial projection feature tensor that perfectly matches the detection channel and angle is found, and the weight value is bound to this spatial position to ensure each weight value accurately corresponds to a spatial position in the initial projection feature tensor. During the matching process, the coordinates of each spatial position are strictly checked to avoid misalignment. After all weight values have been matched, these weight values that are precisely aligned with the spatial positions of the initial projection feature tensor are combined according to their corresponding spatial positions to form spatially aligned weight data.
[0057] Step 302 involves performing element-wise multiplication of the spatial alignment weight data with the corresponding values in the initial projection feature tensor to obtain weighted modulation feature data. Specifically, this includes: clarifying the spatial positional correspondence between the spatial alignment weight data and the initial projection feature tensor, ensuring a one-to-one correspondence between each spatial position; that is, the weight value corresponding to a certain detection channel and angle in the spatial alignment weight data corresponds to the feature value corresponding to the same detection channel and angle in the initial projection feature tensor; performing element-wise multiplication: for each spatial position, extracting the spatial alignment weight value corresponding to that position, and then extracting the corresponding feature value in the initial projection feature tensor; multiplying the extracted weight value and feature value to obtain the weighted modulation value corresponding to that spatial position; specifically, the weighted modulation value equals the spatial alignment weight value multiplied by the feature value at the corresponding position in the initial projection feature tensor; performing the above multiplication operation separately for each spatial position; after all spatial positions have completed the multiplication operation, combining all weighted modulation values according to their corresponding spatial positions to form the weighted modulation feature data.
[0058] Step 303: Perform preset dynamic range boundary truncation processing on the weighted modulation feature data to obtain amplitude-limited enhancement feature data; perform linear sequence reconstruction on the amplitude-limited enhancement feature data to obtain the enhanced projection data stream. Specifically, this includes: determining the preset dynamic range boundary, which is set based on the conventional numerical range of low-dose CT image features, and combining the numerical distribution characteristics of the weighted modulation feature data to set reasonable upper and lower limits. The specific expression for the preset dynamic range boundary is: ,in This represents the preset lower limit value. This represents the preset upper limit value. The value is the mean of all values in the weighted modulation feature data minus 1.5 times the standard deviation of the weighted modulation feature data. The value is the mean of all values in the weighted modulation feature data plus 1.5 times the standard deviation of the weighted modulation feature data. The specific expressions are as follows: ,in This represents the mean of all values in the weighted modulation feature data. This represents the standard deviation of all values in the weighted modulation feature data. The purpose of setting the preset dynamic range boundary is to eliminate outliers and minima caused by noise interference in the weighted modulation feature data, avoiding interference from abnormal values in the subsequent reconstruction process and ensuring the stability of the feature data. After determining the preset dynamic range boundary, each value in the weighted modulation feature data is judged and truncated one by one. If a weighted modulation feature value is greater than the preset upper limit value... Then the value will be truncated to the preset upper limit. If a certain weighted modulation feature value is less than a preset lower limit value Then the value will be truncated to the preset lower limit. If a certain weighted modulation feature value is at a preset upper limit value and preset lower limit value If the value is between these ranges, it remains unchanged and is not truncated. The above judgment and truncation operations are performed individually for each value, ensuring that no value is missed and guaranteeing that all weighted modulation feature values are within the preset dynamic range boundaries. After all the truncation and processing are completed, the data is combined to form amplitude-limited enhancement feature data. This data can effectively remove abnormal interference values and retain valid feature values, further enhancing the reliability of the feature data. The amplitude-limited enhancement feature data is then subjected to linear sequence recombination processing. The recombination method is to arrange all the values in the amplitude-limited enhancement feature data linearly according to the order of the detection channels and the detection angles. Specifically, the arrangement order is as follows: first, arrange the feature values corresponding to all detection angles of the first detection channel, then arrange the feature values corresponding to all detection angles of the second detection channel, and so on, until the feature values corresponding to all detection channels and all detection angles are arranged. During the arrangement process, the spatial position information corresponding to each value is strictly maintained to ensure that each value in the linear sequence can correspond to a specific detection channel and detection angle. After recombination, the enhanced projection data stream is obtained.
[0059] This embodiment, through precise spatial location matching, enables attention weights to be applied precisely to the corresponding positions of the initial projection feature tensor, thereby enhancing the edge features of small lesions and weakening background and noise interference. Furthermore, by dynamically truncating and eliminating abnormal values, it ensures the reliability and effectiveness of the enhanced projection data stream, solving the problem of weak effective features and excessive interference in low-dose CT images.
[0060] In a preferred embodiment of the present invention, step 4, substituting the enhanced projection data stream into the algebraic iterative reconstruction model to obtain a preliminary three-dimensional voxel matrix, may include:
[0061] Step 401: Analyze the scanning geometric trajectory parameters carried in the enhanced projection data stream, extract the set of angle coordinates, and obtain the projection angle parameter sequence; construct a three-dimensional voxel data mesh with all zero assignments based on the spatial imaging area covered by the projection angle parameter sequence to obtain the initial iterative voxel matrix. Specifically, this includes: performing full analysis on the enhanced projection data stream obtained in step 303 using a data analysis algorithm. The specific calculation process of the data analysis algorithm is as follows: first, read the binary header file of the enhanced projection data stream, and extract the total length L of the data stream and the length of the feature value segment. Geometric parameter segment length ,satisfy The reading range of the feature value segment is determined to be [0, The reading range of the geometric parameter segment is []. ,L [1]; Next, the geometric parameter segment data is unpacked, and the binary data of each parameter is read sequentially according to the preset parameter storage order (probe channel coordinates, scanning angle, penetration path coordinates, scanning range, scanning layer thickness, pixel spacing). The binary data is then converted into decimal values. The conversion formula for coordinate parameters is: Decimal value = Binary data × , The coordinate transformation coefficient is set to 0.1 mm to match the accuracy of CT scans. The angle parameter conversion formula is: Decimal angle = Binary data. The angle conversion coefficient m is set to 0.01° to ensure that the parsed parameter values are consistent with the actual scanning parameters. The enhanced projection data stream adopts a binary streaming storage format. The first half stores the feature values after weighted modulation and amplitude limiting processing, while the second half embeds and stores the scanning geometric trajectory parameters synchronously acquired during the low-dose spiral CT scan. These parameters are key imaging information recorded in real time by the scanning device, including the coordinates of the detector channel arrangement, the scanning angle corresponding to each detector channel, the coordinates of the complete X-ray penetration path from the emitter to the detector, the three-dimensional spatial range of the lung scanning area, the scanning slice thickness, and the pixel spacing, as well as all imaging-related geometric information. This information is the core foundation for subsequent construction of the three-dimensional voxel data grid, accurate execution of X-ray penetration integration processing, and ensuring the accuracy of iterative reconstruction. None of these can be omitted. During the parsing process, the feature value segment and the geometric parameter segment of the data stream are first split using the above data parsing algorithm. Then, the geometric parameter segment is parsed in detail to extract the set of angle coordinates. The set of angle coordinates contains all the detector angle values during the low-dose spiral CT scan, covering all angle information from the scan start angle to the scan end angle, without any angle omissions.
[0062] After extraction, these detection angle values are sequentially arranged according to the actual scanning sequence of low-dose spiral CT, forming a projection angle parameter sequence. The arrangement order is completely consistent with the sequence of detection angle changes during the scanning process, that is, arranged sequentially from the first detection angle to the last detection angle according to the angular rotation direction of the scanning equipment. This ensures that the projection angle parameter sequence can completely and accurately reflect the angular changes during the scanning process, providing a precise angular basis for subsequent angle-by-angle X-ray penetration integration processing. After obtaining the projection angle parameter sequence, the size and range of the three-dimensional voxel data grid are further determined based on the spatial imaging area covered by the sequence, where the spatial imaging area is the human body. The complete scanning area of the lungs is determined by the scanning range information in the scanning geometry trajectory parameters. Specifically, it is bounded by the spatial coordinates of the upper and lower edges, left and right edges, and anterior and posterior edges of the lungs, ensuring that the grid can completely cover the entire lung scanning area without missing any lung tissue. The size of the three-dimensional voxel data grid is precisely set according to the actual size of the lung scanning area and the clinical imaging accuracy requirements. The side length of the voxels is combined with the minimum diameter of small lesions to ensure that the size of each voxel can accurately capture the detailed features of small lesions such as ground-glass nodules in the lungs. This avoids the loss of lesion details due to excessively large voxel sizes, or the data redundancy and increased iterative computation burden due to excessively small voxel sizes. After determining the grid size and range, a three-dimensional voxel data grid with all zero values is constructed using a three-dimensional array initialization algorithm. The specific calculation process of the three-dimensional array initialization algorithm is as follows: The three dimensions of the three-dimensional voxel data grid are defined, corresponding to the left-right direction (X-axis), front-back direction (Y-axis), and up-down direction (Z-axis) of the lung scan area, respectively. These three dimensions correspond to the rows, columns, and layers of the three-dimensional array, where the X-axis is the row dimension, the Y-axis is the column dimension, and the Z-axis is the layer dimension. The number of voxels in each dimension is calculated based on the boundary coordinates of the lung scan area, with the number of voxels in the X-axis direction being... Number of voxels in the Y-axis direction Number of voxels in the Z-axis direction ,in , These are the X-axis coordinates of the left and right edges of the lung, respectively. , These are the Y-axis coordinates of the anterior and posterior edges of the lung, respectively. , , where are the Z-axis coordinates of the upper and lower edges of the lung, and 'a' is the voxel side length. The function is used for rounding down, ensuring the number of voxels is an integer and the mesh completely covers the lung region; the 3D array V[ is initialized. ][ ][ The row index range of the array is... The column index range is The range of the layer index is The values of all elements in the array (i.e., each voxel) are set to 0, without assigning any initial feature values, ensuring that the initial iterative voxel matrix is in a state without initial bias; after construction, the initial iterative voxel matrix is obtained, which is... OK, List, The three-dimensional matrix of the layer has three dimensions corresponding to the three spatial dimensions of the lung scanning area: X-axis (left and right), Y-axis (front and back), and Z-axis (up and down). Each voxel corresponds to a tiny spatial region within the lung scanning area, and its spatial position corresponds one-to-one with the spatial position of the actual lung tissue.
[0063] Step 402: Perform ray penetration integration processing on the initial iterative voxel matrix along the current angle coordinate in the projection angle parameter sequence. Calculate the sum of the products of each voxel value and its corresponding penetration length along the penetration path to obtain the forward projection data. Specifically, this includes: selecting the first angle coordinate from the projection angle parameter sequence as the current angle coordinate; subsequently, after each complete iteration update, selecting the next angle coordinate as the current angle coordinate according to the order of the projection angle parameter sequence until all angle coordinates have been selected and the full-angle ray penetration integration processing is completed; after determining the current angle coordinate, based on the geometric information corresponding to that angle coordinate in the scanning geometric trajectory parameters, determine the penetration path of each X-ray at that angle coordinate using a path calculation algorithm. The penetration path is the complete path of the X-ray starting from the emission end, penetrating the human lung scanning area, and finally reaching the corresponding detection channel. Each penetration path and... A detection channel is precisely matched and corresponds one-to-one with an actual projection value in the enhanced projection data stream, ensuring that the forward projection values calculated subsequently can accurately match the actual projection values. The ray penetration integration processing flow is initiated. The specific calculation process is as follows: based on the current angle coordinates and the penetration path, the penetration path is divided into several tiny line segments. The length of each tiny line segment is controlled to be one-tenth of the voxel side length to ensure calculation accuracy. The number of segments is determined by the total length of the penetration path and the length of each tiny line segment, i.e., the total length of the penetration path is divided by the length of each tiny line segment, and rounded up to obtain the minimum number of tiny line segments required. The penetration path is then evenly divided into this number of tiny line segments. The start and end coordinates of each tiny line segment are calculated by uniform distribution, i.e., the start coordinates of each tiny line segment are obtained by adding the product of the tiny line segment length and the line segment number to the start coordinates of the penetration path, and the end coordinates are obtained by adding the start coordinates to the tiny line segment length.
[0064] The voxel to which each tiny line segment belongs is determined using a spatial coordinate judgment algorithm. First, the spatial coordinates of the midpoint of each tiny line segment are calculated. The midpoint coordinates are the average of the coordinates of the starting and ending points of the tiny line segment. Then, based on the starting coordinates of the 3D voxel mesh and the voxel side length, the voxel mesh position to which the midpoint coordinates belong is determined. The method is to calculate the difference between the midpoint coordinates and the mesh starting coordinates, divide it by the voxel side length, and round down to obtain the voxel mesh index in that direction. The combination of the indices in the three directions gives the mesh coordinates of the voxel to which the tiny line segment belongs. The penetration length corresponding to each voxel is then calculated. The penetration length corresponding to each voxel is the sum of the lengths of all the tiny line segments within that voxel. By counting the tiny line segments contained in each voxel one by one and summing the lengths of these tiny line segments, the actual penetration length of that voxel on the penetration path is obtained. The current value of each voxel is extracted, and the current value of each voxel is multiplied by the penetration length corresponding to that voxel to obtain the product value of each voxel. Then, the product values of all voxels are summed to obtain the forward projection value corresponding to the penetration path, thus completing the ray penetration integral calculation for a single penetration path. For each penetration path, the ray penetration integral is calculated one by one. Take all voxels along the path, determine the specific position of each voxel on the penetration path, and then, based on the voxel side length and the angle between the penetration path and the voxel, calculate the corresponding penetration length of the voxel through geometric calculation. This is calculated by multiplying the cosine of the angle between the penetration path and the voxel by the voxel side length, yielding the actual length of the penetration path within the voxel. This length represents the actual physical distance the X-ray travels through the voxel. During the calculation, accurately calculate the influence of the angle on the penetration length to avoid calculation errors. For each penetration path, perform numerical multiplication and accumulation operations, first extracting the value of each voxel on the penetration path... The current value is then multiplied by the penetration length corresponding to that voxel. Each voxel undergoes this multiplication operation individually, without batch merging, to ensure the accuracy of each product value. After multiplying all voxels, the product values of all voxels along the penetration path are accumulated one by one, following the order in which the voxels appear on the path, ensuring no voxel's product value is missed. The sum is the forward calculated projection value corresponding to that penetration path, reflecting the theoretical projection value of the corresponding detection channel at the current angular coordinates. Each penetration path undergoes this multiplication and accumulation operation independently to ensure that the calculation results of different penetration paths are independent and do not interfere with each other. The forward calculated projection values corresponding to all penetration paths are combined according to the order of the detection channels to form the forward calculated projection data.
[0065] Step 403: Extract the actual projection values corresponding to the current angle coordinates from the enhanced projection data stream. Perform point-by-point subtraction between the forward-calculated projection data and the actual projection values to obtain the projection residual data sequence. Specifically, this includes: based on the current angle coordinate values, extracting the actual projection values corresponding to all detection channels at that angle coordinate from the enhanced projection data stream using a coordinate matching algorithm. The specific calculation process of the coordinate matching algorithm is as follows: extract the scanning geometric trajectory parameters embedded in the enhanced projection data stream, clarify the numbering range of the detection channels, the angle coordinate range corresponding to each detection channel, and the linear recombination rule of the enhanced projection data stream, i.e., linearly arranged according to the detection channel order and the order of all angle coordinates corresponding to each channel; calculate the index position of the current angle coordinate in the angle sequence, and determine the index position based on the arrangement rule of the enhanced projection data stream. The current angular coordinates correspond to the starting and ending offsets of all detection channel values. The starting offset equals the current angular coordinate index multiplied by the total number of detection channels, and the ending offset equals the starting offset plus the total number of detection channels, ensuring accurate positioning of the numerical range corresponding to the current angular coordinates. Following the order of the detection channel numbers, each value within the range is extracted one by one. During extraction, the correspondence between the detection channel number and the value is strictly verified; that is, detection channel number 1 corresponds to the first value in the range, detection channel number 2 corresponds to the second value, and so on, until all values corresponding to all detection channels have been extracted. Each extracted value is then bound to its corresponding detection channel number and the current angular coordinates, forming a three-dimensional correspondence between the detection channel, angular coordinates, and actual projected value, thus completing coordinate matching and value extraction.
[0066] The actual projection value is the feature value in the enhanced projection data stream obtained in step 303, after weighted modulation, amplitude limiting, and linear sequence recombination, that precisely matches the current angle coordinates and the corresponding detection channel. It is also the actual scan projection value corresponding to the detection channel at that angle coordinate. During the extraction process, the number of each detection channel and the value of the current angle coordinate are strictly checked to ensure that the extracted actual projection value corresponds precisely to the current angle coordinates and detection channel, without any angle misalignment, channel misalignment, or value extraction error, thus avoiding distortion of subsequent residual calculations due to extraction errors. After extraction, the forward calculated projection data and the extracted actual projection value are subtracted point by point. The point-by-point subtraction strictly follows the one-to-one correspondence principle. For each detection channel at the current angle coordinate, the forward calculated projection value corresponding to that detection channel in the forward calculated projection data is extracted first, and then... The actual projection value corresponding to the detection channel is extracted from the enhanced projection data stream. The forward-calculated projection value is subtracted from the actual projection value to obtain the projection residual value corresponding to the detection channel. Specifically, the projection residual value is equal to the actual projection value minus the forward-calculated projection value. If the actual projection value is greater than the forward-calculated projection value, the residual value is positive; if the actual projection value is less than the forward-calculated projection value, the residual value is negative; if the two are equal, the residual value is zero. The magnitude of the residual value directly reflects the degree of deviation between the theoretical projection value and the actual projection value. The above subtraction operation is performed individually for each detection channel, without omitting any detection channel or cross-subtraction between different detection channels, ensuring the accuracy of each projection residual value. The projection residual values corresponding to all detection channels are arranged one by one in the order of the detection channel numbers to form a projection residual data sequence.
[0067] Step 404: Residual values in the projection residual data sequence are back-mapped to the spatial grid positions corresponding to the initial iterative voxel matrix according to the penetration path of the current angle coordinates. Weight allocation calculation is performed to obtain the voxel update increment matrix. Specifically, this includes: clarifying the correspondence between each penetration path under the current angle coordinates and the spatial grid positions in the initial iterative voxel matrix; confirming the spatial grid coordinates of all voxels corresponding to each penetration path through a spatial coordinate back-mapping algorithm. The specific calculation process of the spatial coordinate back-mapping algorithm is as follows: extracting the complete spatial coordinate information of each penetration path under the current angle coordinates, including the starting coordinates, ending coordinates, and endpoint coordinates of all small line segments obtained by discretization in step 402, forming a complete set of penetration path coordinates; calling the starting coordinates of the three-dimensional voxel data grid in step 401. The coordinates, voxel side lengths, and voxel quantity parameters in each dimension are used to define the spatial range of the 3D voxel mesh and the spatial coordinate interval of each voxel. For each penetration path, each coordinate point in its coordinate set is selected one by one. The voxel mesh coordinates to which the coordinate point belongs are calculated using spatial coordinate judgment logic. The calculation method is consistent with the voxel positioning logic in step 402, that is, the difference between the coordinate point and the starting coordinate of the mesh is calculated, then divided by the voxel side length, and rounded down to obtain the mesh index in three directions. The combination is the mesh coordinate of the voxel to which the coordinate point belongs. The voxel mesh coordinates to which all coordinate points corresponding to each penetration path belong are deduplicated to remove duplicate mesh coordinates, so as to obtain the mesh coordinates of all unique voxels traversed by the penetration path, forming a one-to-one correspondence between the penetration path and the voxel mesh coordinates, and completing the spatial coordinate reverse mapping. Each penetration path corresponds to a series of continuous voxels in the initial iterative voxel matrix, which is completely consistent with and one-to-one with the penetration path processed by the ray penetration integration in step 402, ensuring that the residual value can be accurately reverse mapped to the corresponding voxel position without any mapping misalignment. Each residual value in the projection residual data sequence is back-mapped to the spatial grid position of the initial iterative voxel matrix according to the corresponding penetration path. That is, the residual value corresponding to each detection channel is back-mapped to all voxels on the penetration path corresponding to that detection channel. One residual value corresponds to multiple voxels, ensuring that the deviation of the residual value can be evenly distributed to all voxels on the penetration path, realizing the accurate transmission of deviation. After the mapping is completed, a weight allocation calculation is performed on each voxel. The core basis of the weight allocation is the penetration length of the voxel on the penetration path. The principle of the longer the penetration length is followed, the greater the weight is, because the longer the penetration length, the greater the influence of the voxel on the projection value, and the greater the corresponding residual allocation weight should be. Conversely, the smaller the weight is, the more the weight allocation can fit the actual physical process of X-ray penetration of voxels, making the residual correction more targeted.The specific calculation process is as follows: First, for each penetration path, the total penetration length of all voxels along that path is calculated. The penetration length of each voxel along that path is then divided by the total penetration length of that path to obtain the weight allocation coefficient for each voxel. The weight allocation coefficient ranges from 0 to 1, and the sum of the weight allocation coefficients for all voxels is 1 to ensure the rationality of the weight allocation. Next, the residual value corresponding to the penetration path is multiplied by the weight allocation coefficient for each voxel to obtain the update increment value for that voxel. Specifically, the update increment value equals the residual value multiplied by the weight allocation coefficient. When the residual value is positive, the update increment value is positive, and the voxel value will increase; when the residual value is negative, the update increment value is negative, and the voxel value will decrease. The magnitude of the increment value reflects the degree to which the voxel needs correction. The above weight allocation calculation is performed individually for each voxel, avoiding batch allocation or average allocation, to ensure that the update increment value for each voxel is accurate and reasonable, and can specifically correct the deviation between the voxel and the actual lung tissue. The update increment values of all voxels are combined one by one according to their corresponding spatial grid positions to form a voxel update increment matrix. The dimension of this matrix is exactly the same as the dimension of the initial iteration voxel matrix, and each voxel position corresponds to a unique update increment value.
[0068] Step 405 involves performing point-by-point addition between the voxel values in the initial iterative voxel matrix and the corresponding increment values in the voxel update increment matrix to obtain the current iterative voxel matrix. Specifically, this includes confirming the spatial correspondence between the initial iterative voxel matrix and the voxel update increment matrix. Both matrices have identical three-dimensional dimensions, and each voxel's spatial grid coordinates correspond one-to-one. That is, the voxel corresponding to a specific spatial coordinate in the initial iterative voxel matrix precisely matches the voxel corresponding to the same spatial coordinate in the voxel update increment matrix, ensuring accurate point-by-point addition without misalignment. The point-by-point addition process uses independent voxel-by-voxel calculations for each voxel. For each voxel in the voxel matrix, its current value is first extracted using its coordinate index. This value is either the zero value after initialization of the initial iteration voxel matrix or the value updated in the previous iteration. Then, using the same coordinate index, the update increment value at the corresponding position in the voxel update increment matrix is extracted. The current value and the update increment value are added together to obtain the updated value of the voxel. Specifically, the updated value of the voxel equals the value of the voxel in the initial iteration voxel matrix plus the increment value of the voxel in the voxel update increment matrix. The addition operation is strictly followed to ensure calculation accuracy and prevent numerical overflow or calculation errors. The above addition operation is performed on each voxel individually, without omitting any voxel. At the same time, it is strictly guaranteed that the current value of each voxel and the increment value at the corresponding position are precisely matched, without cross-addition or misaligned addition, to avoid voxel value distortion due to calculation errors, which would affect the reconstruction effect of subsequent iterations. After all voxels have been added, the current iteration voxel matrix is formed.
[0069] Step 406: Calculate the global numerical change rate between the current iteration voxel matrix and the previous iteration reference matrix. Update the current iteration voxel matrix to the next iteration reference matrix. Repeat the ray penetration integration process to point-by-point addition until the global numerical change rate is lower than the convergence threshold to obtain the preliminary three-dimensional voxel matrix. Specifically, this includes: defining the previous iteration reference matrix, where in the first iteration, the previous iteration reference matrix is the initial iteration voxel matrix constructed in step 401, with all voxel values being zero; starting from the second iteration, the previous iteration reference matrix is the current iteration voxel matrix obtained in the previous iteration process, i.e., the voxel matrix corrected in the previous iteration, which serves as the comparison reference for this iteration; initiate the calculation process of the global numerical change rate, strictly following the quantitative comparison principle. For each voxel in the current iteration voxel matrix, extract the current value of the voxel through a one-to-one coordinate correspondence, and then extract the reference value of the voxel corresponding to the same coordinate in the previous iteration reference matrix. The numerical difference is obtained by subtracting two values. The absolute value of each voxel's difference is then taken to eliminate the influence of negative differences and ensure their validity. This yields the absolute difference for each voxel. The absolute differences of all voxels in the current iteration voxel matrix are then summed, with each voxel's absolute difference counted sequentially to obtain a global sum. This sum reflects the overall deviation between the current iteration voxel matrix and the previous iteration's baseline matrix. Finally, the sum of the absolute values of all voxel values in the previous iteration's baseline matrix is calculated, and then divided by the global sum of differences to obtain the global rate of change. Specifically, the global rate of change equals the global sum of differences divided by the sum of the absolute values of voxel values in the previous iteration's baseline matrix. This value ranges from 0 to 1; a smaller value indicates a smaller deviation between the current iteration voxel matrix and the previous iteration's baseline matrix, and the iteration effect is closer to convergence.
[0070] After the calculation is completed, the current iteration voxel matrix is updated to the reference matrix for the next iteration, overwriting the reference matrix of the previous iteration, providing a new comparison reference for the next iteration. All operations from steps 402 to 405 are repeated, that is, the next current angle coordinate is selected again from the projection angle parameter sequence, and the following operations are performed in sequence: ray penetration integration processing, forward calculation of projection data and subtraction of actual projection values, reverse mapping of residual values and weight allocation calculation, and point-by-point addition of voxel values to obtain a new current iteration voxel matrix. Then, the global numerical change rate between the new current iteration voxel matrix and the updated reference matrix of the previous iteration is calculated. This process is repeated iteratively, with each iteration specifically correcting the deviation of the voxel matrix and gradually improving the matching degree between the voxel matrix and the actual lung tissue, until the calculated global numerical change rate is lower than the preset convergence threshold. The convergence threshold is set strictly in line with the clinical accuracy requirements of low-dose CT image reconstruction and is set in combination with the accuracy requirements for identifying small lung lesions to ensure that the voxel matrix after iteration can accurately reflect the actual lung imaging situation. At this point, all iteration operations are stopped, and the voxel matrix obtained in the last iteration is determined as the preliminary three-dimensional voxel matrix.
[0071] In this embodiment, an initial voxel matrix is constructed by analyzing the scanning geometric trajectory parameters. Then, by combining the ray penetration integral, residual calculation, reverse mapping weight allocation and iterative update, the deviation between the voxel matrix and the actual imaging is gradually corrected, so that the preliminary three-dimensional voxel matrix can accurately reflect the three-dimensional structure of the lung and clearly present the details of tiny lesions.
[0072] In a preferred embodiment of the present invention, step 5, extracting the three-dimensional local neighborhood features of the preliminary three-dimensional voxel matrix and constructing a three-dimensional volumetric attention mask; using the three-dimensional volumetric attention mask to perform anisotropic diffusion smoothing and detail compensation on the preliminary three-dimensional voxel matrix to obtain a three-dimensional reconstructed image, may include:
[0073] Step 501: Analyze the mesh topology of the preliminary three-dimensional voxel matrix, extract the spatial position index of each voxel node, and obtain the three-dimensional spatial coordinate sequence. Specifically, this includes: reading the mesh configuration information of the preliminary three-dimensional voxel matrix, clarifying the three-dimensional mesh topology of the matrix, and the mesh topology specifically including the number of voxels in each dimension, the voxel side length, and the starting spatial coordinates of the three-dimensional voxel mesh. The number of voxels in each dimension is consistent with the initial iterative voxel matrix constructed in step 401, corresponding to the left-right, front-back, and up-down directions of the lung scanning area, respectively. The voxel side length is also consistent with the voxel side length set in step 401, ensuring that the spatial correspondence between the mesh topology and the actual lung scanning area remains unchanged. During the analysis of the mesh topology, the number and side length of voxels in each dimension were checked one by one to ensure that the mesh could completely cover the human lung scanning area without any missing mesh or exceeding the scanning area. Simultaneously, the continuity of the mesh topology was ensured to avoid mesh breaks or misalignments, providing a stable mesh foundation for subsequent extraction of voxel node spatial position indices. After analysis, the spatial position indices of each voxel node were extracted. Each voxel node's spatial position index consists of three index values, corresponding to the row index i, column index j, and layer index k of the initial 3D voxel matrix. Row index i corresponds to the left and right sides of the lung. In the direction extraction, column index j corresponds to the anterior-posterior direction of the lung, and layer index k corresponds to the superior-inferior direction of the lung. During extraction, the spatial position index of each voxel node is extracted one by one in the following order: row index i from 0 to the number of voxels in the left-right direction minus 1, column index j from 0 to the number of voxels in the anterior-posterior direction minus 1, and layer index k from 0 to the number of voxels in the superior-inferior direction minus 1. This ensures that no voxel node is missed and that there are no duplicate extractions. After extracting the spatial position index of each voxel node, the actual spatial coordinates corresponding to each spatial position index are calculated by combining the starting spatial coordinates of the 3D voxel mesh and the voxel side length. Let the starting coordinates of the left-right direction be... The starting coordinates in the forward and backward directions are The starting coordinates in the up and down directions are If the side length of the voxel is *a*, then the actual coordinates in each direction are calculated as follows: Actual coordinates x in the left and right directions = +i×a, actual coordinates y in the forward and backward directions = +j×a, actual vertical coordinates z= +k×a; binds the spatial position index of each voxel node to its corresponding actual spatial coordinates (x,y,z), and arranges them one by one according to the spatial distribution order of the voxel nodes to form a three-dimensional spatial coordinate sequence.
[0074] Step 502: Based on the three-dimensional spatial coordinate sequence, a cube neighborhood is defined centered on the voxel node corresponding to each spatial location index. The set of voxel values contained within the cube neighborhood is extracted to obtain the three-dimensional local neighborhood features. Specifically, this includes: defining the size setting standard for the cube neighborhood; the side length of the cube neighborhood is set according to the voxel side length and the size of the small lung lesion. Let the voxel side length be a, and the cube neighborhood side length be n×a, where n is an odd number. This ensures that the neighborhood centered on the target voxel node is symmetrical, improving the accuracy of local feature extraction, while avoiding the neighbor range being too large, causing interference from irrelevant voxel values, or too small, failing to fully capture local structural features. After determining the cube neighborhood size, for each voxel node in the three-dimensional spatial coordinate sequence, the actual spatial coordinates of that voxel node ( , , Using a cube as the center, define its neighborhood by extending it in six directions: left, right, front, back, up, and down, with a distance of [missing information]. The boundary coordinates of the cube's spatial range are calculated as follows: lower limit in the left and right directions. = Upper limit = Lower limit in the forward and backward directions = Upper limit = Lower limit in the vertical direction = Upper limit = The upper and lower limits of the six directions together constitute the spatial range of the cube's neighborhood. After defining the cube's neighborhood, all voxel nodes contained within that neighborhood are extracted. The extraction process involves checking the actual spatial coordinates (x, y, z) of each voxel node in the preliminary 3D voxel matrix one by one to determine whether they satisfy the following conditions. ≤x≤ , ≤y≤ , ≤z≤ If all conditions are met, the voxel node belongs to the neighborhood of the cube and is selected. If any direction does not meet the condition, it is not extracted. After extraction, the voxel values corresponding to all voxel nodes within the neighborhood of the cube are collected to form a voxel value set. Each voxel value set is bound to the corresponding central voxel node. For all central voxel nodes, the above operation is performed one by one. After the neighborhood features of all voxel nodes are extracted, all voxel value sets are combined one by one according to the spatial position order of their corresponding central voxel nodes to obtain the three-dimensional local neighborhood features.
[0075] Step 503: Based on the three-dimensional local neighborhood features, calculate the local variance and spatial gradient aggregation value within the neighborhood of each cube to obtain the set of local structural distribution parameters. Specifically, this includes: for the set of voxel values corresponding to each cube neighborhood in the three-dimensional local neighborhood features, calculate the local variance within that neighborhood. The local variance reflects the dispersion of the voxel values within that neighborhood. The greater the dispersion, the more obvious the changes in lung tissue structure within that neighborhood, and the possible presence of micro-lesions; the smaller the dispersion, the more uniform the lung tissue structure within that neighborhood, and the more likely it is normal tissue. The specific calculation process for the local variance is as follows: first, calculate the average value of all voxel values within the cube neighborhood. The calculation method is to add all voxel values within the neighborhood to obtain the sum of voxel values, and then divide the sum of voxel values by... The average voxel value is obtained by using the total number of voxel nodes in the neighborhood. Then, the difference between each voxel value and the average value is calculated, and each difference is squared to obtain the squared difference of each voxel value. All squared differences are summed to obtain the sum of squared differences. Finally, the sum of squared differences is divided by the total number of voxel nodes in the neighborhood to obtain the local variance of the cube neighborhood. The local variance ranges from 0 to 1000, which is consistent with the distribution characteristics of voxel values in low-dose CT. After calculating the local variance of each cube neighborhood, the spatial gradient aggregation value within each cube neighborhood is calculated. The spatial gradient aggregation value is used to reflect the spatial variation trend of voxel values in the neighborhood, which can capture the abrupt changes in voxel values in the neighborhood, and thus identify the edge of lung tissue and the outline of small lesions. The specific calculation process for the spatial gradient aggregation value is as follows: For each voxel node in the neighborhood of the cube, calculate the numerical difference between the voxel node and its neighboring voxel nodes. Neighboring voxel nodes include those in the six directions (up, down, left, right, front, and back). Calculate the numerical difference between the target voxel and each of its neighboring voxels, and take the absolute value of each difference to obtain the absolute difference in the six directions. Then, sum the absolute differences in the six directions corresponding to the voxel node to obtain the spatial gradient value of that voxel node. Sum the spatial gradient values of all voxel nodes in the neighborhood of the cube to obtain the sum of the spatial gradients of the neighborhood. Finally, divide the sum of the spatial gradients by the total number of voxel nodes in the neighborhood to obtain the spatial gradient of the cube's neighborhood. The local variance and spatial gradient aggregation values range from 0 to 500. Each cube neighborhood is calculated individually, avoiding batch calculations to ensure the accuracy of each parameter calculation and prevent errors. After calculation, the local variance and spatial gradient aggregation values for each cube neighborhood are bound together and associated with the spatial location index of the corresponding central voxel node, ensuring that each central voxel node has a unique local variance and spatial gradient aggregation value. All local variance and spatial gradient aggregation values corresponding to the central voxel nodes are then combined sequentially according to their spatial location to form a set of local structural distribution parameters.
[0076] Step 504: Input each parameter in the local structure distribution parameter set into a preset nonlinear mapping function to perform numerical transformation, obtaining initial volume weight values. Combine all initial volume weight values in spatial order to obtain initial volume weight distribution data. Perform three-dimensional spatial neighborhood mean filtering on the initial volume weight distribution data to obtain smooth volume weight distribution data. Specifically, this includes: setting the input and output ranges of the nonlinear mapping function. The input range corresponds to the conventional value range of local variance and spatial gradient aggregation values, and the output range is set between 0 and 1 to ensure that the initial volume weight values can be used as weight coefficients. The preset nonlinear mapping function adopts an improved form of the Sigmoid function, and its mathematical expression is: ,in The input parameters are local variance or spatial gradient aggregation value. `c` is the gain coefficient, used to adjust the rate of change of the function, ranging from 0.001 to 0.01. `b` is the threshold parameter, used to adjust the offset of the function, ranging from the mean of the local variance or spatial gradient aggregation value, ensuring that the function output matches the weighting requirements. The construction of the nonlinear mapping function must meet two core requirements: first, the larger the local variance and the larger the spatial gradient aggregation value, the larger the initial volume weight value of the output, because these areas are likely to be small lung lesions or tissue edges, requiring higher weights and subsequent optimization; second, the smaller the local variance and the smaller the spatial gradient aggregation value, the smaller the initial volume weight value of the output, because these areas are likely to be normal, uniform tissue, and the weights can be relatively low. After construction, the nonlinear mapping function is validated by selecting local variance and spatial gradient aggregation values of different ranges as input to the function, checking whether the initial volume weight value of the output meets the preset requirements. If not, the function parameters `c` and `b` are adjusted until the output result meets the requirements, completing the construction of the preset nonlinear mapping function.
[0077] Each local variance and its corresponding aggregated spatial gradient value in the set of local structure distribution parameters are input into the preset nonlinear mapping function to perform numerical transformation. The specific transformation process is as follows: Let the local variance be... ², The corresponding weight components are obtained by inputting the nonlinear mapping function. =f( ²), weighted components The value ranges from 0 to 0.5; let the spatial gradient aggregation value be g, and input the nonlinear mapping function to obtain the corresponding weight components. =f(g), weighted components The value ranges from 0 to 0.5; the two weight components are added together to obtain the initial volume weight value corresponding to the central voxel node. The specific calculation process is as follows: initial volume weight value w = + The initial volume weight value w ranges from 0 to 1, ensuring that the initial volume weight value can comprehensively reflect the dispersion and spatial variation trend of the local structure.
[0078] Each central voxel node undergoes the aforementioned numerical transformation operation individually to obtain its corresponding initial volume weight value. All initial volume weight values are then combined sequentially according to their spatial positions within the corresponding central voxel nodes to form the initial volume weight distribution data. This data has the same dimension as the initial 3D voxel matrix, with each voxel position corresponding to an initial volume weight value, accurately reflecting the importance of each voxel node. Since the initial volume weight distribution data may contain abrupt changes in weight values due to fluctuations in local structural parameters, affecting the stability of subsequent mask construction, a 3D spatial neighborhood mean filtering process is required. The specific filtering process involves setting the size of the 3D spatial neighborhood, which is consistent with the size of the cube neighborhood in step 502, ensuring that the filtering process is effective. Based on local structural features, for each voxel node in the initial volume weight distribution data, a three-dimensional spatial neighborhood is defined with that voxel node as the center, and the initial volume weight values corresponding to all voxel nodes within that neighborhood are extracted. Then, the average value of all initial volume weight values within that neighborhood is calculated by adding all initial volume weight values within the neighborhood to obtain the weight sum, and then dividing the weight sum by the total number of voxel nodes within the neighborhood to obtain the mean value. Finally, this mean value is used as the filtered weight value of the central voxel node, replacing the original initial volume weight value. The smoothed volume weight value still ranges from 0 to 1. The above filtering operation is performed on all voxel nodes individually. After filtering all voxel nodes, all filtered weight values are combined in order of their corresponding spatial positions to obtain the smoothed volume weight distribution data.
[0079] Step 505: Extract the global maximum and global minimum values from the smoothed volume weight distribution data. Map the initial volume weight values to a preset confidence interval according to the range ratio to obtain a 3D volume attention mask. Specifically, this includes: traversing all weight values in the smoothed volume weight distribution data, comparing the size of each weight value one by one, and selecting the largest weight value as the global maximum value of the smoothed volume weight distribution data. The global maximum value ranges from 0.8 to 1. At the same time, select the smallest weight value as the global minimum value of the smoothed volume weight distribution data. The global minimum value ranges from 0 to 0.2. Ensure that no value is overlooked during the extraction process. To avoid errors in extracting the global maximum or minimum value, any weight value must be omitted, ensuring the accuracy of subsequent scaling. After extraction, the weight range of the smoothed volume weight distribution data is calculated. The weight range is calculated by subtracting the global minimum value from the global maximum value, and the weight range range is from 0.6 to 1. A preset confidence interval is set according to the weight requirements of the 3D volume attention mask, usually between 0 and 1, to ensure that the mapped weight values can be used as mask weights, conforming to the weight value specifications in the attention mechanism, and meeting the weight requirements for subsequent anisotropic diffusion smoothing and detail compensation. Each initial volume weight value in the initial volume weight distribution data is mapped to a pre-set confidence interval according to the weight range. Specifically, the mapping process involves calculating the difference between each initial volume weight value and the global minimum value to obtain a single weight difference; then dividing this single weight difference by the weight range to obtain a scaling factor; multiplying this scaling factor by the length of the pre-set confidence interval (where the interval is 1 if the interval is between 0 and 1), the scaling factor is directly used as the mapped weight value. If the initial volume weight value is less than the global minimum value, the mapped weight value is the minimum value of the pre-set confidence interval, 0; if the initial volume weight value is less than the global minimum value, the mapped weight value is the minimum value of the pre-set confidence interval, 0. If the value is greater than the global maximum, the mapped weight value is taken as the maximum value of the preset confidence interval, which is 1. Ensure that the mapped weight values are all within the preset confidence interval of 0 to 1, while preserving the relative differences of the initial volume weight values. Perform the above proportional mapping operation on each initial volume weight value separately to obtain the mapped volume weight value. Combine all the mapped volume weight values one by one according to the spatial position of their corresponding voxel nodes to form a three-dimensional volume attention mask. The dimension of this mask is completely consistent with the dimension of the initial three-dimensional voxel matrix and the initial volume weight distribution data. Each voxel position corresponds to a unique weight value.
[0080] Step 506: Spatial position matching and alignment of the weight values in the 3D volumetric attention mask with the corresponding voxel values in the preliminary 3D voxel matrix to obtain mask-aligned voxel data; based on the mask-aligned voxel data, extracting the difference between the voxel values corresponding to each spatial position and the values of adjacent voxels in a preset neighborhood to obtain the voxel spatial gradient field. Specifically, this includes: clarifying the spatial position correspondence between the 3D volumetric attention mask and the preliminary 3D voxel matrix, ensuring that the dimensions of the two are completely consistent, and that each voxel position corresponds one-to-one, that is, the weight value at a certain spatial position in the 3D volumetric attention mask corresponds to the voxel value at the same spatial position in the preliminary 3D voxel matrix, ensuring the accuracy of the matching and alignment, and keeping the weight values of the 3D volumetric attention mask within the range of 0 to 1. A spatial alignment operation is performed. For each voxel location, the weight value corresponding to that location in the 3D volumetric attention mask is extracted, and the corresponding voxel value in the preliminary 3D voxel matrix is also extracted. The extracted weight values and voxel values are then bound together. During the binding process, the spatial positions are strictly checked to avoid misalignment and ensure that each voxel value accurately matches its corresponding weight value, preventing situations where one voxel value corresponds to multiple weight values or vice versa. After all voxel locations are aligned, all voxel values bound with weight values are combined one by one according to their corresponding spatial positions to obtain mask-aligned voxel data. This data contains both the voxel value information from the preliminary 3D voxel matrix and the weight information from the 3D volumetric attention mask, providing fused data support for subsequent gradient calculation and diffusion smoothing. Based on the mask-aligned voxel data, the difference between the voxel values at each spatial location and the values of adjacent voxels within a preset neighborhood is extracted to obtain the voxel spatial gradient field. A preset neighborhood is set, with its size consistent with the size of the cube neighborhood in step 502, ensuring that the neighborhood range fits the local structural features. The preset neighborhood includes adjacent voxel nodes in six directions (up, down, left, right, front, and back) for each voxel node, comprehensively capturing the spatial changes in voxel values. For each voxel node in the mask-aligned voxel data, the voxel value corresponding to that node is extracted, and then extracted one by one... Extract the voxel values corresponding to the six adjacent voxel nodes within its preset neighborhood. Calculate the difference in value between the target voxel node and each adjacent voxel node. Specifically, the difference in value is equal to the voxel value of the target voxel node minus the voxel values of the adjacent voxel nodes. Calculate the difference in value for each direction separately to obtain the difference in value for all six directions. A positive difference indicates that the target voxel value is greater than the adjacent voxel values; a negative difference indicates that the target voxel value is less than the adjacent voxel values; and a zero difference indicates that the target voxel value is equal to the adjacent voxel values.Each voxel node performs the above-mentioned variation difference extraction operation individually to obtain the variation difference in six directions corresponding to that voxel node. The variation differences in six directions corresponding to all voxel nodes are combined one by one in order of their corresponding spatial positions to form the voxel spatial gradient field.
[0081] Step 507: Based on the voxel spatial gradient field and the weight values in the 3D volumetric attention mask, calculate the diffusion conduction coefficient corresponding to each voxel position to obtain a spatially non-uniform diffusion coefficient distribution; based on the spatially non-uniform diffusion coefficient distribution, perform neighborhood weighted iterative update processing on the voxel values in the mask-aligned voxel data to obtain a smoothed 3D voxel matrix. Specifically, this includes: clarifying the calculation logic of the diffusion conduction coefficient, which is used to control the intensity of anisotropic diffusion smoothing. In regions with larger gradients and larger weights, the diffusion conduction coefficient is smaller, and the diffusion smoothing intensity is weaker, ensuring the preservation of lesion edges and details; in regions with smaller gradients and smaller weights, the diffusion conduction coefficient is larger, and the diffusion smoothing intensity is stronger, achieving noise suppression of normal tissue areas, meeting the needs of noise reduction and edge preservation in low-dose CT lung reconstruction.
[0082] When calculating the diffusion conduction coefficient corresponding to each voxel position, for each voxel node, the difference in variation in six directions in the voxel space gradient field is extracted. This difference in variation is essentially the partial derivative of the voxel value in each direction in three-dimensional space. The partial derivatives in the six directions are expressed as follows: (Left and right directions) (Reverse left and right directions) (Front-back direction) (Reverse front and back direction) (up and down direction), (Inverted up and down directions), where I is the voxel value, and x, y, and z correspond to the left-right, front-back, and up-down directions in 3D space, respectively; calculate the absolute value of each partial derivative to obtain the absolute value difference in the six directions, and then add the six absolute value differences to obtain the sum of the absolute values of the gradient of that voxel node, calculated as G = Where G is the sum of absolute gradients; simultaneously, the weight values of the voxel node in the 3D volumetric attention mask are extracted. , The value ranges from 0 to 1. The sum of the absolute values of the gradients G and the weight values are calculated together. Multiplying them together yields the gradient weight product, which is calculated as P = G × Where P is the gradient weight product; a fixed smoothing coefficient is set. , The value ranges from 0.1 to 0.3, and is set according to noise reduction requirements and edge preservation accuracy to ensure that the diffusion smoothing effect meets expectations; finally, the smoothing coefficient is... Dividing by the product of gradient weights P and the sum of 1, we obtain the diffusion conduction coefficient D corresponding to the voxel node. The specific calculation expression is D = The diffusion conduction coefficient D ranges from 0.07 to 0.3, ensuring it is between 0 and 1. The larger the gradient weight product P, the smaller the diffusion conduction coefficient D, and vice versa.
[0083] Each voxel node performs the aforementioned diffusion conduction coefficient calculation operation independently, obtaining the diffusion conduction coefficient corresponding to each voxel position. All diffusion conduction coefficients are then combined sequentially according to their corresponding spatial positions to obtain a spatially non-uniform diffusion coefficient distribution. This distribution enables differentiated diffusion control for different lung regions, providing accurate coefficient support for subsequent neighborhood-weighted iterative updates. Based on the spatially non-uniform diffusion coefficient distribution, neighborhood-weighted iterative update processing is performed on the voxel values in the mask-aligned voxel data. First, the number of iterations is set according to the smoothing effect requirements, ensuring effective noise suppression without losing lung tissue details. For each voxel node in the mask-aligned voxel data, in each iteration, a preset neighborhood range is defined centered on that voxel node. This preset neighborhood is consistent with the preset neighborhood set in step 506. The voxel values and diffusion conduction coefficients corresponding to all adjacent voxel nodes within this neighborhood range are extracted.
[0084] The weighted value of each adjacent voxel node is calculated as follows: the weighted value equals the voxel value of the adjacent voxel node multiplied by the diffusion conduction coefficient of that adjacent voxel node, and then multiplied by the corresponding 3D volumetric attention mask weight value. The 3D volumetric attention mask weight value ranges from 0 to 1. The weighted values of all adjacent voxel nodes within the neighborhood are summed to obtain the total weighted value. Simultaneously, the total weight coefficient of all adjacent voxel nodes within the neighborhood is calculated. The total weight coefficient is equal to the sum of the diffusion conduction coefficient of each adjacent voxel node multiplied by its corresponding weight value. The total weight coefficient ranges from 0 to the volumetric attention mask weight value within the neighborhood. The total number of voxel nodes; if the sum of weight coefficients is not zero, the sum of weighted values is divided by the sum of weight coefficients to obtain the iterative update value of the voxel node; if the sum of weight coefficients is zero, the iterative update value of the voxel node remains unchanged and is still the current voxel value; the current voxel value of the voxel node is replaced with the iterative update value to complete one iterative update; all voxel nodes perform the above neighborhood weighted iterative update operation individually, and after completing one iteration, they enter the next iteration until the preset number of iterative updates is reached; after the iterative update is completed, the final iterative update values of all voxel nodes are combined in order of their corresponding spatial positions to obtain the smoothed three-dimensional voxel matrix.
[0085] Step 508: Based on the smoothed 3D voxel matrix, extract the gradient jump trajectory of the 2D cross-section layer by layer along the 3D spatial axis to identify the local structural missing regions and obtain the local structural missing contour set; based on the local structural missing contour set, extract the internal spatial coordinate nodes enclosed by the local structural missing contour set, and perform edge-constrained voxel numerical interpolation filling processing on the internal spatial coordinate nodes to obtain the locally connected repair 3D voxel matrix. Specifically, this includes: clearly defining the division of the 3D spatial axis, which are the left-right axis, the front-back axis, and the up-down axis, and extracting the 2D cross-section layer by layer along the three axes to ensure that all spatial regions of the smoothed 3D voxel matrix are fully covered without missing any 2D cross-section, which fits the complexity of the lung's 3D structure and ensures that the local structural missing regions can be fully identified. When extracting two-dimensional cross sections layer by layer along each axis, the two-dimensional cross section data of the corresponding layer is extracted one by one according to the index order of the axis. Each two-dimensional cross section data contains the voxel values of all voxel nodes in that layer. During the extraction process, the spatial position correspondence of voxel nodes is strictly maintained to avoid cross section misalignment or value loss. After extracting each two-dimensional cross section data, the gradient jump trajectory within the cross section is calculated. The gradient jump trajectory is used to identify regions where voxel values change abruptly within the cross section. These regions are usually the edges of lung tissue. If the edges are interrupted or discontinuous, it corresponds to a region with missing local structure. The extraction process of the gradient jump trajectory is as follows: for each voxel node in the two-dimensional cross section, the numerical difference between the voxel node and its adjacent voxel nodes is calculated. Adjacent voxel nodes include adjacent voxels in both horizontal and vertical directions. The numerical difference in each direction is calculated, and the absolute value is taken. If the absolute value is greater than a preset gradient threshold, it indicates that the voxel node is at a gradient jump position, and the spatial coordinates of the voxel node are recorded. After traversing all voxel nodes in the two-dimensional cross section, the spatial coordinates of all gradient jump positions are connected in sequence to form the gradient jump trajectory of the two-dimensional cross section. The preset gradient threshold is set according to the numerical distribution characteristics of the smoothed three-dimensional voxel matrix, with a value range of 5 to 20, to ensure that the real tissue edge can be accurately identified and to avoid misjudging small numerical changes caused by noise as gradient jumps.After extracting the gradient transition trajectories of all two-dimensional sections, a coherence analysis is performed on each gradient transition trajectory. If a gradient transition trajectory is interrupted, broken, or discontinuous, it indicates a local structural deficiency in that region. The gradient transition trajectory segment corresponding to the interrupted region is extracted as the local structural deficiency contour. Simultaneously, by combining the gradient transition trajectories of the two-dimensional sections along the three axes, the local structural deficiency contours of different sections are spatially integrated to form a local structural deficiency contour in three-dimensional space, ensuring accurate location of the spatial extent of the local structural deficiency. After all local structural deficiency contours are integrated, they are combined according to their spatial position order to obtain a set of local structural deficiency contours, which clearly reflects... The process involves smoothing all regions with missing local structures in the 3D voxel matrix to provide accurate contour information for subsequent repair processing. Based on the set of missing local structure contours, the internal spatial coordinate nodes enclosed by each missing local structure contour are extracted. The extraction process involves determining the spatial region enclosed by each missing local structure contour and checking the spatial coordinates of voxel nodes in the smoothed 3D voxel matrix one by one. If the spatial coordinates of a voxel node are within the internal region enclosed by the missing local structure contour, then the voxel node is an internal spatial coordinate node and is extracted. If it is outside or on the contour, it is not extracted, ensuring that all extracted internal spatial coordinate nodes are voxel nodes within the missing local structure region.
[0086] After extraction, edge-constrained voxel numerical interpolation and filling are performed on the internal spatial coordinate nodes. The core of interpolation and filling is to constrain the interpolated values of the internal spatial coordinate nodes based on the edge voxel values of the missing local structural contours, ensuring that the filled voxel values are consistent with the voxel values of the surrounding normal tissue and conform to the natural structure of the lung tissue. Specifically, for each internal spatial coordinate node, the edge voxel nodes on the missing local structural contours around the node are extracted, and the voxel values of the edge voxel nodes are extracted. The average value of these edge voxel values is calculated by adding all the edge voxel values to obtain the sum of the edge values, and then dividing the sum of the edge values by the total number of edge voxel nodes to obtain the average edge value. At the same time, the spatial distance between the internal spatial coordinate node and each edge voxel node is calculated. Take the reciprocal of each spatial distance to obtain a distance weight, with the weight ranging from 0.1 to 1. The closer the distance, the larger the weight; the farther the distance, the smaller the weight. Multiply each edge voxel value by its corresponding distance weight to obtain a weighted edge value. Sum all weighted edge values to obtain a total weighted edge value. Then sum all distance weights to obtain a total distance weight. Divide the total weighted edge value by the total distance weight to obtain the interpolated fill value for that internal spatial coordinate node. Replace the original voxel value of that internal spatial coordinate node with the interpolated fill value to complete the filling process for that node. Perform the above interpolation filling process on all internal spatial coordinate nodes individually. After filling, combine all the repaired voxel values according to their corresponding spatial positions to obtain a locally connected repair 3D voxel matrix.
[0087] Step 509: Calculate the numerical difference between the locally connected 3D voxel matrix and the smoothed 3D voxel matrix to obtain high-frequency detail residual data. Based on a 3D volumetric attention mask, perform amplitude enhancement processing on the components corresponding to high-weight values in the high-frequency detail residual data to obtain compensated detail data. Specifically, this includes: clarifying the spatial correspondence between the locally connected 3D voxel matrix and the smoothed 3D voxel matrix, ensuring that their dimensions are completely identical and each voxel position corresponds one-to-one to guarantee the accuracy of the numerical difference calculation. Calculate the numerical difference between the two matrices. For each voxel position, extract the voxel value corresponding to that position in the locally connected 3D voxel matrix, and then extract the voxel value corresponding to that position in the smoothed 3D voxel matrix. Subtract the voxel value of the smoothed 3D voxel matrix from the voxel value of the locally connected 3D voxel matrix to obtain the numerical difference corresponding to that voxel position. Specifically, the numerical difference is equal to the locally connected voxel value minus the smoothed voxel value. A positive numerical difference indicates that the value at that voxel location has increased after repair, supplementing missing details; a negative numerical difference indicates that the value at that voxel location has decreased after repair, correcting redundant values; a zero numerical difference indicates that there is no structural loss at that voxel location, and the value remains unchanged. The above numerical difference calculation operation is performed individually on all voxel locations to obtain the numerical difference corresponding to each voxel location. All numerical differences are combined one by one according to their corresponding spatial positions to obtain high-frequency detail residual data. This data contains high-frequency detail information supplemented during the local structural repair process. This information mainly consists of details of small lung lesions and subtle features of tissue edges, which is crucial for improving the accuracy of 3D reconstructed images. Based on a 3D volumetric attention mask, amplitude enhancement processing is performed on the components corresponding to high-weight values in the high-frequency detail residual data. A high-weight value threshold is set, which is set according to the weight distribution characteristics of the 3D volumetric attention mask, with a value range of 0.6 to 0.8, ensuring accurate screening of high-weight regions, i.e., regions with a high probability of small lesions and tissue edges. The filtering process involves iterating through each weight value in the 3D volumetric attention mask. The weight value ranges from 0 to 1. If a weight value is greater than a preset high weight value threshold, the area corresponding to that weight value is a high weight area, and the spatial location of that area is recorded. If the weight value is less than or equal to the preset high weight value threshold, the area is a low weight area, and no amplitude enhancement processing is performed.After screening, amplitude enhancement processing is performed on the numerical differences corresponding to high-weight regions in the high-frequency detail residual data. Specifically, a fixed enhancement coefficient is set, ranging from 1.2 to 1.8, based on the detail compensation requirements, to ensure that the enhanced details are clearly visible without artifacts caused by over-enhancement. The numerical difference corresponding to each voxel position in the high-weight region is multiplied by the enhancement coefficient to obtain the enhanced numerical difference. Specifically, the enhanced numerical difference equals the original numerical difference multiplied by the enhancement coefficient. The numerical differences in low-weight regions remain unchanged and are not enhanced to avoid interfering with the details of normal tissue areas. After amplitude enhancement processing is completed for all high-weight regions, the enhanced numerical differences are combined with the unenhanced numerical differences according to their corresponding spatial positions to obtain the compensated detail data.
[0088] Step 510 involves performing point-by-point overlay processing on the compensated detail data and the voxel values in the smoothed 3D voxel matrix to obtain a detail-enhanced voxel matrix. Spatial coordinate mapping and grayscale standardization are then performed on the detail-enhanced voxel matrix to obtain the 3D reconstructed image. Specifically, this includes: clarifying the spatial correspondence between the compensated detail data and the smoothed 3D voxel matrix, ensuring their dimensions are completely identical and each voxel position corresponds one-to-one; that is, the numerical difference at a certain spatial position in the compensated detail data corresponds to the voxel value at the same spatial position in the smoothed 3D voxel matrix, ensuring the accuracy of point-by-point overlay. Point-by-point overlay processing is performed. For each voxel position, the numerical difference corresponding to that position in the compensated detail data is extracted, and then the corresponding voxel value in the smoothed 3D voxel matrix is extracted. The smoothed voxel value is added to the difference in compensated detail values to obtain the detail-enhanced voxel value corresponding to that voxel position. Specifically, the detail-enhanced voxel value equals the smoothed voxel value plus the difference in compensated detail values.
[0089] The above stacking operation is performed individually at each voxel location, without batch stacking, to ensure stacking accuracy and avoid computational errors. Simultaneously, the difference between the smoothed voxel value and the compensated detail value at each voxel location is strictly matched to prevent cross-stacking. After stacking at all voxel locations, all detail-enhancing voxel values are combined sequentially according to their corresponding spatial locations to obtain a detail-enhancing voxel matrix. This matrix retains the noise reduction effect of the smoothed 3D voxel matrix while supplementing enhanced high-frequency detail information, clearly presenting the normal tissue structure and minute lesion details of the lungs, thus improving the imaging accuracy of the 3D voxel matrix. Spatial coordinate system mapping is performed on the detail-enhancing voxel matrix. The purpose of spatial coordinate system mapping is to map the grid coordinate system of the detail-enhancing voxel matrix to the clinically commonly used medical imaging coordinate system to ensure the accuracy of the 3D reconstructed image. To adapt to medical imaging diagnostic equipment, the specific mapping process involves defining the origin, axis direction, and scaling ratio of the medical imaging coordinate system. The scaling ratio ranges from 0.01 to 0.1 millimeters per pixel. The mesh coordinates of each voxel node in the detail enhancement voxel matrix are converted into actual spatial coordinates in the medical imaging coordinate system. The conversion process involves the medical imaging coordinates in each direction being equal to the mesh coordinates multiplied by the voxel side length, plus the origin offset of the corresponding direction in the medical imaging coordinate system. This ensures that the mapped coordinates accurately correspond to the actual anatomical location of the human lung, meeting the needs of clinical diagnosis. After spatial coordinate system mapping, grayscale range standardization is performed on the detail enhancement voxel matrix. The purpose of grayscale range standardization is to adjust the grayscale range of voxel values to the standard grayscale range of clinical medical images, ensuring clear grayscale contrast between different regions. Specifically, this process involves extracting the maximum and minimum grayscale values of all voxel values in the detail enhancement voxel matrix, calculating the grayscale range (equal to the maximum grayscale value minus the minimum grayscale value), setting a standard grayscale range based on the conventional grayscale range of clinical medical images (typically 0 to 255), and mapping each voxel value to the standard grayscale range according to its grayscale range. The specific mapping process is as follows: Subtracting the minimum gray value from the voxel value yields the gray-level difference. Dividing this difference by the gray-level range gives a scaling factor, which ranges from 0 to 1. Multiplying this factor by the length of the standard gray-level interval and adding the minimum value of that interval gives the standardized gray-level value. If the voxel value is less than the minimum gray value, the standardized gray-level value is the minimum of the standard gray-level interval. If the voxel value is greater than the maximum gray value, the standardized gray-level value is the maximum of the standard gray-level interval. All voxel values undergo this gray-level interval standardization process individually. After processing, all standardized gray-level values are combined in order of their corresponding spatial locations to obtain the 3D reconstructed image.
[0090] This embodiment can accurately capture the structural features of local lung tissue, especially the local features of small lesions, by extracting the three-dimensional local neighborhood features of the preliminary three-dimensional voxel matrix, which helps to achieve differentiated treatment of different regions of the lung.
[0091] like Figure 2 As shown, embodiments of the present invention also provide a CT image three-dimensional reconstruction system based on an attention mechanism, comprising:
[0092] The parsing module is used to collect sparse view sinusoidal data, parse and normalize the sparse view sinusoidal data, and construct an initial projection feature tensor containing high-frequency edge gradient information and low-frequency background baseline information.
[0093] The computation module is used to input the initial projected feature tensor into a pre-constructed multi-head self-attention encoding network to obtain the attention query feature tensor and the attention key feature tensor; and to calculate the dot product similarity between the transposes of the attention query feature tensor and the attention key feature tensor to obtain a two-dimensional spatial attention weight map.
[0094] The processing module is used to perform element-wise weighted modulation processing on the initial projection feature tensor based on the two-dimensional spatial attention weight map to obtain the enhanced projection data stream;
[0095] The reconstruction module is used to substitute the enhanced projection data stream into the algebraic iterative reconstruction model to obtain a preliminary three-dimensional voxel matrix.
[0096] The output module is used to extract the three-dimensional local neighborhood features of the preliminary three-dimensional voxel matrix and construct a three-dimensional volumetric attention mask. The three-dimensional volumetric attention mask is used to perform anisotropic diffusion smoothing and detail compensation on the preliminary three-dimensional voxel matrix to obtain a three-dimensional reconstructed image.
[0097] It should be noted that this system is a system corresponding to the above method. All implementation methods in the above method embodiments are applicable to this embodiment and can achieve the same technical effect.
[0098] Embodiments of the present invention also provide a computing device, including: a processor and a memory storing a computer program, wherein the computer program, when executed by the processor, performs the method described above. All implementations in the above method embodiments are applicable to this embodiment and can achieve the same technical effects.
[0099] The above description represents the preferred embodiments of the present invention. It should be noted that those skilled in the art can make various improvements and modifications without departing from the principles of the present invention, and these improvements and modifications should also be considered within the scope of protection of the present invention.
Claims
1. A 3D reconstruction method for CT images based on an attention mechanism, characterized in that, The method includes: Step 1: Collect sparse view sinusoidal data, analyze and standardize the sparse view sinusoidal data, and construct an initial projection feature tensor containing high-frequency edge gradient information and low-frequency background baseline information. Step 2: Input the initial projection feature tensor into the pre-constructed multi-head self-attention encoding network to obtain the attention query feature tensor and the attention key feature tensor; calculate the dot product similarity between the transposes of the attention query feature tensor and the attention key feature tensor to obtain the two-dimensional attention weight map. Step 3: Based on the two-dimensional spatial attention weight map, the initial projection feature tensor is subjected to element-wise weighted modulation processing to obtain the enhanced projection data stream; Step 4: Substitute the enhanced projection data stream into the algebraic iterative reconstruction model to obtain the preliminary three-dimensional voxel matrix; Step 5: Extract the 3D local neighborhood features of the preliminary 3D voxel matrix and construct a 3D volumetric attention mask; use the 3D volumetric attention mask to perform anisotropic diffusion smoothing and detail compensation on the preliminary 3D voxel matrix to obtain a 3D reconstructed image, including: analyzing the mesh topology of the preliminary 3D voxel matrix, extracting the spatial position index of each voxel node, and obtaining a 3D spatial coordinate sequence; based on the 3D spatial coordinate sequence, delineate the cube neighborhood range with the voxel node corresponding to each spatial position index as the center, extract the set of voxel values contained within the cube neighborhood range, and obtain the 3D local neighborhood features; based on the 3D local neighborhood features, calculate the local variance and spatial gradient aggregation value within the cube neighborhood range to obtain the local structure analysis. The local structure distribution parameter set is configured; each parameter in the local structure distribution parameter set is input into a preset nonlinear mapping function to perform numerical transformation, obtaining initial volume weight values. All initial volume weight values are combined in spatial order to obtain initial volume weight distribution data. Three-dimensional spatial neighborhood mean filtering is performed on the initial volume weight distribution data to obtain smooth volume weight distribution data. The global maximum and global minimum values are extracted from the smooth volume weight distribution data, and the initial volume weight values are mapped to a preset confidence interval according to the range ratio to obtain a three-dimensional volume attention mask. The weight values in the three-dimensional volume attention mask are spatially matched and aligned with the voxel values at corresponding spatial positions in the initial three-dimensional voxel matrix to obtain the mask. Align voxel data; based on the mask-aligned voxel data, extract the difference between the voxel values at each spatial location and the values of adjacent voxels in a preset neighborhood to obtain the voxel spatial gradient field; based on the voxel spatial gradient field and the weight values in the 3D volume attention mask, calculate the diffusion conduction coefficient corresponding to each voxel location to obtain the spatially non-uniform diffusion coefficient distribution; based on the spatially non-uniform diffusion coefficient distribution, perform neighborhood-weighted iterative update processing on the voxel values in the mask-aligned voxel data to obtain a smoothed 3D voxel matrix; based on the smoothed 3D voxel matrix, extract the 2D cross-sectional gradient jump trajectory layer by layer along the 3D spatial axis to identify local structural missing regions and obtain a set of local structural missing contours; based on the set of local structural missing contours... The process involves extracting internal spatial coordinate nodes enclosed by a set of missing local structural contours, performing edge-constrained voxel numerical interpolation on these nodes to obtain a locally connected 3D voxel matrix. The numerical difference between the locally connected 3D voxel matrix and the smoothed 3D voxel matrix is calculated to obtain high-frequency detail residual data. Based on a 3D volumetric attention mask, amplitude enhancement processing is performed on the components corresponding to high-weight values in the high-frequency detail residual data to obtain compensated detail data. The compensated detail data and the voxel values in the smoothed 3D voxel matrix are then superimposed point-by-point to obtain a detail-enhanced voxel matrix. Finally, spatial coordinate system mapping and grayscale normalization are performed on the detail-enhanced voxel matrix to obtain the 3D reconstructed image.
2. The CT image three-dimensional reconstruction method based on attention mechanism according to claim 1, characterized in that, Sparse viewpoint sinusoidal image data is collected, and the data is analyzed and normalized to construct an initial projection feature tensor containing high-frequency edge gradient information and low-frequency background baseline information, including: Acquire sparse view sinusoidal plot data, perform logarithmic domain transformation and channel consistency compensation processing on the sparse view sinusoidal plot data to obtain linearized corrected projection data; Based on linearized corrected projection data, the neighborhood numerical difference of each sampling point in the preset detection direction is extracted, and the gradient magnitude aggregation processing of the neighborhood numerical difference is performed to obtain high-frequency edge gradient information. The low-frequency background baseline information is obtained by subtracting the high-frequency edge gradient information at the corresponding position from the values of each sampling point in the linearized corrected projection data. Based on high-frequency edge gradient information, the global maximum and global minimum values are extracted, and the values of each sampling point in the high-frequency edge gradient information are mapped to a preset standard interval according to the range ratio to obtain standardized high-frequency edge gradient information. Based on the low-frequency background baseline information, the global maximum and global minimum values are extracted, and the values of each sampling point in the low-frequency background baseline information are mapped to a preset standard interval according to the range ratio to obtain the standardized low-frequency background baseline information. The standardized high-frequency edge gradient information and the standardized low-frequency background baseline information are concatenated and stitched together to obtain the initial projection feature tensor.
3. The CT image three-dimensional reconstruction method based on attention mechanism according to claim 2, characterized in that, The initial projected feature tensor is input into a pre-constructed multi-head self-attention encoding network to obtain the attention query feature tensor and the attention key feature tensor, including: The initial projection feature tensor is input into a pre-constructed multi-head self-attention encoding network, and the initial projection feature tensor is linearly mapped using the network's preset query weight parameters to obtain the initial query feature data. The initial query feature data is divided into multiple independent query feature fragments by performing a grouping and segmentation operation; the multiple independent query feature fragments are then spliced and recombined according to a preset data arrangement order to obtain the attention query feature tensor. The initial projection feature tensor is input into a pre-constructed multi-head self-attention encoding network. The initial projection feature tensor is linearly mapped using the network's preset key weight parameters to obtain initial key feature data. The initial key feature data is then divided into multiple independent key feature segments by a grouping and segmentation operation. Multiple independent key feature fragments are spliced and recombined according to a preset data arrangement order to obtain the attention key feature tensor.
4. The CT image three-dimensional reconstruction method based on attention mechanism according to claim 3, characterized in that, Calculate the dot product similarity between the transpose of the attention query feature tensor and the attention key feature tensor to obtain a two-dimensional attention weight map, including: The spatial dimension and feature dimension of the attention key feature tensor are interchanged to obtain the transposed key feature tensor. The attention query feature tensor and the transposed key feature tensor are multiplied and summed at their corresponding spatial locations to obtain the original similarity data. The original similarity data is divided by a preset ratio coefficient to obtain scaled similarity data; the scaled similarity data is then mapped by an exponential function to obtain positive exponential similarity data. The positive exponential similarity data is accumulated position by position along the preset feature dimension to obtain the normalized denominator data corresponding to each spatial position; the value of each spatial position in the positive exponential similarity data is divided by the corresponding normalized denominator data to obtain a two-dimensional spatial attention weight map.
5. The CT image three-dimensional reconstruction method based on attention mechanism according to claim 4, characterized in that, Based on a two-dimensional spatial attention weight map, the initial projection feature tensor is subjected to element-wise weighted modulation processing to obtain an enhanced projection data stream, including: The spatial attention weight map of the two-dimensional space is matched according to the spatial position distribution of the initial projection feature tensor to obtain spatial alignment weight data. The spatially aligned weight data is multiplied element-wise with the corresponding values in the initial projected feature tensor to obtain the weighted modulation feature data. The weighted modulation feature data is subjected to a preset dynamic range boundary truncation process to obtain amplitude-limited enhancement feature data; the amplitude-limited enhancement feature data is then linearly recombined to obtain an enhanced projection data stream.
6. The CT image three-dimensional reconstruction method based on attention mechanism according to claim 5, characterized in that, Substituting the enhanced projection data stream into the algebraic iterative reconstruction model yields a preliminary three-dimensional voxel matrix, including: The scanning geometric trajectory parameters carried in the enhanced projection data stream are analyzed, and the set of angle coordinates is extracted to obtain the projection angle parameter sequence. Based on the spatial imaging area covered by the projection angle parameter sequence, a three-dimensional voxel data grid with all zero assignments is constructed to obtain the initial iterative voxel matrix. The initial iterative voxel matrix is subjected to ray penetration integration along the current angle coordinate in the projection angle parameter sequence. The sum of the products of each voxel value and the corresponding penetration length on the penetration path is calculated to obtain the forward projection data. Extract the actual projection values corresponding to the current angular coordinates from the enhanced projection data stream, and perform point-by-point subtraction between the forward calculated projection data and the actual projection values to obtain the projection residual data sequence. The residual values in the projection residual data sequence are back-mapped to the spatial grid position corresponding to the initial iterative voxel matrix according to the penetration path of the current angle coordinates, and the weight allocation calculation is performed to obtain the voxel update increment matrix. The voxel values in the initial iterative voxel matrix are added point by point to the corresponding increment values in the voxel update increment matrix to obtain the current iterative voxel matrix. Calculate the global numerical rate of change between the current iteration voxel matrix and the reference matrix of the previous iteration, update the current iteration voxel matrix to the reference matrix of the next iteration, and repeat the ray penetration integration process to point-by-point addition process until the global numerical rate of change is lower than the convergence threshold to obtain the preliminary three-dimensional voxel matrix.
7. A CT image three-dimensional reconstruction system based on an attention mechanism, wherein the system implements the method as described in any one of claims 1 to 6, characterized in that, include: The parsing module is used to collect sparse view sinusoidal data, parse and normalize the sparse view sinusoidal data, and construct an initial projection feature tensor containing high-frequency edge gradient information and low-frequency background baseline information. The computation module is used to input the initial projection feature tensor into a pre-constructed multi-head self-attention encoding network to obtain the attention query feature tensor and the attention key feature tensor. Calculate the dot product similarity between the transpose of the attention query feature tensor and the attention key feature tensor to obtain a two-dimensional attention weight map. The processing module is used to perform element-wise weighted modulation processing on the initial projection feature tensor based on the two-dimensional spatial attention weight map to obtain the enhanced projection data stream; The reconstruction module is used to substitute the enhanced projection data stream into the algebraic iterative reconstruction model to obtain a preliminary three-dimensional voxel matrix. The output module is used to extract the three-dimensional local neighborhood features of the preliminary three-dimensional voxel matrix and construct a three-dimensional volumetric attention mask. The three-dimensional volumetric attention mask is used to perform anisotropic diffusion smoothing and detail compensation on the preliminary three-dimensional voxel matrix to obtain a three-dimensional reconstructed image.
8. A computing device, characterized in that, include: One or more processors; A storage device for storing one or more programs that, when executed by one or more processors, cause the one or more processors to implement the method as described in any one of claims 1 to 6.
Citation Information
Patent Citations
Medical image focus segmentation and three-dimensional reconstruction method and system based on artificial intelligence
CN121504944A
Depressed scar three-dimensional shape reconstruction method and system carrying AI deep learning
CN122049298A