Photon counting CT image noise reduction method based on multi-channel WGAN

Abstract

The invention discloses a photon counting CT image noise reduction method based on multi-channel WGAN, and relates to the technical field of medical image processing, and the method comprises the steps: taking high-energy, low-energy and all-round CT images as generator input, extracting multi-scale features through independent convolution branch, carrying out the global average pooling of each channel feature map, and calculating the attention weight between channels to carry out the weighted fusion; a WGAN structure is adopted, and a discriminator is constrained through a Wasserstein distance and a gradient penalty term so as to improve the training stability; meanwhile, the mean square error, the structural similarity index and the perception loss are integrated, and the numerical precision and the visual effect of the noise reduction image are balanced; and by alternately optimizing the generator and the discriminator, performing multi-channel feature processing and up-sampling on the noisy low-dose CT image after training is completed, and outputting a high-quality noise-reduced image.

Description

Technical Field

[0001] The present invention relates to the field of medical image processing technology, and in particular to a photon counting CT image denoising method based on a multi-channel WGAN. Background Art

[0002] In the field of medical image processing, low-dose photon counting CT (PCCT) holds significant clinical diagnostic value due to its low radiation dose and multi-energy spectral imaging advantages. However, low-dose imaging is susceptible to noise, which can degrade image quality. Traditional denoising methods, such as BM3D and Total Variation (TV) filtering, rely on fixed prior assumptions and suffer from issues such as blurred edges and texture loss. Single-channel generative adversarial networks (GANs), while superior in denoising, fail to fully utilize PCCT's multi-energy spectral information (high-energy, low-energy, and full-energy images). Furthermore, conventional GAN ​​models are prone to mode collapse during training, resulting in unstable denoising results.

[0003] Based on the above problems, the present invention proposes a photon counting CT image denoising method based on multi-channel WGAN to solve one or more of the above problems. Summary of the Invention

[0004] The purpose of the present invention is to provide a photon counting CT image denoising method based on multi-channel WGAN to solve the problems raised in the prior art.

[0005] To achieve the above object, the present invention provides the following technical solutions:

[0006] The photon counting CT image denoising method based on multi-channel WGAN includes the following steps:

[0007] S1. Photon counting CT is used to collect high-energy, low-energy, and full-energy images. After noise reduction, correction, and normalization preprocessing, a three-channel input tensor is formed to provide a standardized data basis for subsequent feature extraction.

[0008] S2, extract deep feature maps from the three-channel images through independent convolution branches, and each layer is batch normalized and ReLU activated to obtain multi-scale semantic features;

[0009] S3. Perform global average pooling on the deep feature maps of each channel to generate a global vector, calculate the channel similarity weight through self-attention, and perform weighted fusion to form a fused feature map that focuses on the key features of multiple energy domains;

[0010] S4. Design a multi-layer convolutional neural network as the discriminator, use the Wasserstein distance to measure the distribution difference between the real and generated samples, and constrain it to meet the 1-Lipschitz condition through the gradient penalty term;

[0011] S5. Integrate Wasserstein loss, pixel error, structural similarity and perceptual loss to construct the total loss function of the generator to balance the numerical, structural and semantic consistency of the denoised image;

[0012] S6. Alternately train the generator and discriminator. Fix the generator to optimize the discriminator to improve the sample discrimination ability. Then fix the discriminator to optimize the generator to minimize the total loss. Iterate until the model converges.

[0013] S7. The preprocessed low-dose noisy CT image is input into the trained generator. After multi-channel feature extraction, self-attention fusion and upsampling, a high-quality image with noise suppression and detail preservation is output.

[0014] The step S1 further includes the following contents:

[0015] S1.1: The high-energy image H, low-energy image L, and full-energy image T of the same scanned object are acquired by a photon counting CT device, including the clean image NDCT obtained by normal-dose scanning and the noisy image LDCT obtained by low-dose scanning. Both are two-dimensional grayscale images with a uniform size of h×w, i.e., height×width. Among them, the clean image NDCT obtained by normal-dose scanning refers to the photon counting CT image obtained by standard radiation dose scanning, which is considered as a real image without significant noise after preprocessing; the noisy image LDCT obtained by low-dose scanning refers to the photon counting CT image obtained by low-dose radiation scanning. The image, including noise due to insufficient dose, serves as the noisy input data for the generator. NDCT stands for normal-dose photon-counting CT, and LDCT stands for low-dose CT. Based on the multi-energy spectral imaging characteristics of PCCT, high-energy images reflect the tissue's absorption characteristics of high-energy X-rays, low-energy images correspond to low-energy X-ray absorption characteristics, and full-energy images integrate information from the entire energy range. PCCT is a photon-counting CT device. The acquired images are then preprocessed, including removing background noise, correcting for detector response differences, and normalizing pixel values ​​to the neural network input range, which is [0,1].

[0016] S1.2: The preprocessed high-energy, low-energy, and full-energy images are used as a set of inputs to the generator, and combined with an optional random noise vector Z to form an input tensor with 3 channels; specifically expressed as [H, L, T, Z], where H∈R h×w , L∈R h×w , T∈R h×w, corresponding to the two-dimensional grayscale matrices of high-energy, low-energy, and full-energy images, respectively, h and w are the height and width of the image; the number of input channels of the generator is 3, which does not include channel Z, which is an optional additional input, and then directly receives the tensor as input; the generator will subsequently perform independent feature extraction on each channel, using a 5×5 convolution kernel for the high-energy channel, a 7×7 convolution kernel for the low-energy channel, and a 3×3 convolution kernel for the full-energy channel to ensure that the construction of the input tensor maintains the independence and correspondence of the three channels.

[0017] The step S2 further includes the following contents:

[0018] S2.1: In view of the different characteristic features of high-energy, low-energy and full-energy images, an independent convolutional neural network branch is designed for each channel, namely the high-energy image channel, low-energy image channel and full-energy image channel; for the high-energy image channel, a 5×5 convolution kernel is used, with a total of 8 convolution layers; the number of convolution kernels in the first layer is determined according to the number of input channels, initially set to 3 channels, and the number of convolution kernels in each subsequent layer is gradually increased according to the network design, with 64 convolution kernels in the first layer, and each subsequent layer increases by multiples, with the final number of kernels being 512, which is used to extract mid-scale features and balance noise and detail information; for the low-energy image channel, a 7×7 convolution kernel is used, also containing 8 convolution layers; used to extract generalized noise patterns and global structural features contained in low-energy images; for the full-energy image channel, a 3×3 convolution kernel is used, with 8 convolution layers; the specific operation process of each convolution layer is as follows: after the input image is calculated by the convolution kernel, the features are normalized by the batch normalization layer, and then nonlinear mapping is introduced by the ReLU activation function;

[0019] S2.2: The 8 convolutional layers of each channel extract features from the input image layer by layer. The front convolutional layer captures shallow features, and the back convolutional layer gradually integrates them to form deep semantic features. The 5×5 kernel of the high-energy channel focuses on extracting medium-scale noise-signal mixed features in the middle convolutional layer. The 7×7 kernel of the low-energy channel covers a large area in the shallow layer to capture the noise distribution law. The 3×3 kernel of the full-energy channel focuses on retaining fine textures in the deep layer. Finally, each channel outputs the corresponding deep feature map F. H 、F L 、F T Among them, F H 、F L 、F T They are represented as high-energy, low-energy and full-energy deep feature maps respectively.

[0020] The step S3 further includes the following contents:

[0021] S3.1: For each channel deep feature map F generated in step S2 H 、F L 、F TPerform global average pooling operations respectively; global average pooling operation is referred to as GAP operation later; these feature maps are all three-dimensional tensors with dimension R C×h×w , where C is the number of channels of the feature map, which is determined by the configuration of the generator convolutional layer, and h×w is the spatial size of the feature map, that is, the height × width of the image; the GAP operation compresses the three-dimensional tensor into a one-dimensional global description vector v∈R by performing channel-by-channel mean calculation on the spatial dimension of each feature map. C , the specific calculation is:

[0022] v H =GAP(F H )∈R C , v L =GAP(F L )∈R C , v T =GAP(F T )∈R C ;

[0023] in,

[0024]

[0025] c∈{1,2,...,C} is the channel index, that is, all spatial pixels of each channel are averaged to generate the global statistical features of the channel;

[0026] S3.2: Based on the one-dimensional global description vector v obtained in sub-step 3.1 H ,v L ,v T , calculate the attention weights between channels through the self-attention mechanism; transform the one-dimensional global description vector vH,vL,vT∈R C Through the fully connected layer mapping to the query matrix Qi and the key matrix Kj, i, j∈{H,L,T}, calculate the similarity matrix A∈R 3×3 , where A i,j =Q i ·K j , characterizes the pairwise correlation between the three channels; the query matrix Qi and the key matrix Kj, i, j∈{H,L,T} correspond to high energy, low energy, and full energy channels respectively; the dimension after mapping remains R C ; Then calculate the similarity between channel i and channel j by matrix dot product to obtain the similarity matrix A∈R 3×3 , where A i,j =Q i ·K j ; Then perform Softmax processing on the similarity matrix to generate normalized attention weights. The formula is as follows:

[0027]

[0028] Among them, α H , α L , α T are the attention weights of high-energy, low-energy, and full-energy channels, respectively, and α H +α L +α T =1, indicating the contribution of each channel in feature fusion;

[0029] S3.3: Use the normalized attention weight α calculated in substep 3.2 H , α L , α T , for each channel deep feature map F H 、F L 、F T Perform weighted fusion, and the dimensions of each feature map are R C×h×w , the fusion calculation formula is as follows:

[0030] F fused =α H ·F H +α L ·F L +α T ·F T ;

[0031] Among them, F fused ∈R C×h×w The pixel value of each channel of the fused feature map is calculated by the weighted sum of the corresponding positions of the original feature map; the contribution of each channel feature is dynamically adjusted through the attention weight, and the generator will focus on the key features of different energy domains: the generalized noise pattern extracted by the 7×7 convolution kernel in the low-energy channel, the edge details extracted by the 5×5 kernel in the high-energy channel, and the texture information extracted by the 3×3 kernel in the all-energy channel.

[0032] The step S4 further includes the following contents:

[0033] S4.1: Design the network architecture of the discriminator. Construct it as a 6-layer convolutional neural network. Each layer uses a 3×3 convolution kernel. The configuration of each layer is as follows:

[0034] For the convolution layer configuration, the first two layers are 64 convolution kernels with stride 1 and 64 convolution kernels with stride 2, respectively. The next two layers are 128 convolution kernels with stride 1 and 128 convolution kernels with stride 2, and the next two layers are 256 convolution kernels with stride 1 and 256 convolution kernels with stride 2. Each layer is followed by a LeakyReLU activation function. The convolution operation with stride 2 reduces the spatial size of the input image by half according to the preset range. At the same time, the number of channels is doubled layer by layer, and the low-level pixel features are gradually extracted to the high-level semantic features.

[0035] For the fully connected layer configuration, after 6 convolution layers, a 1024-dimensional fully connected layer and a 1-dimensional output layer are sequentially connected to output the score of the input sample. This structure does not use batch normalization and meets the Wasserstein GAN's 1-Lipschitz constraint on the discriminator.

[0036] S4.2: The discriminator receives a real clean image sample x and a denoised image sample G(z) output by the generator. Image sample x comes from a normal-dose photon-counting CT scan, a noise-free image acquired and preprocessed in step S1.1; denoised image sample G(z) comes from a noisy low-dose CT image input to the generator in step S1.2 and processed in steps S2-S3. The approximate Wasserstein distance between the two is then calculated. The first two terms of the loss function are as follows:

[0037] L(D)=-E x～Preal [D(x)]+E z～Pnoise [D(G(z))];

[0038] Among them, L(D) is the loss function of the discriminator, which is used to measure the approximate Wasserstein distance between the true distribution and the generated distribution. E is the mathematical expectation, which represents the average operation of the sample distribution. D(x) is the score of the discriminator on the true sample x. D(G(z)) is the score of the discriminator on the generated sample G(z). Preal is the distribution of the NDCT true image, and Pnoise is the distribution of the LDCT noisy image. By maximizing L(D), the discriminator is forced to improve the score of the true sample and reduce the score of the generated sample, thereby guiding the generator to adjust the parameters so that the distribution of G(z) approaches the true distribution p. real ;

[0039] Then constrain the discriminator to meet the 1-Lipschitz condition and uniformly sample intermediate samples along the line connecting the real sample and the generated sample Calculate the gradient penalty term as follows:

[0040]

[0041] Where λ is a constant weighting parameter; Represents the distribution of the middle sample Expected operation of ; Represents the discriminator output for the intermediate sample The gradient vector reflects the slope of the discriminator in the sample space; Represents the L2 norm of the gradient vector, which measures the length of the gradient; Represents the core of the penalty term, forcing the gradient norm to be close to 1 so that the discriminator satisfies the 1-Lipschitz constraint;

[0042] S4.3: Perform adversarial training, alternately optimizing the discriminator and generator. Specifically: input a batch of real samples x and generated samples G(z), calculate the total loss L(D), which includes the Wasserstein distance term and the gradient penalty term, and update the discriminator parameters through backpropagation. The goal is to maximize L(D) and improve its ability to distinguish between real and generated samples; then rely on the gradient penalty term to constrain the discriminator gradient norm.

[0043] The step S5 further includes the following contents:

[0044] S5.1: Pixel-wise error L between the denoised image G(z) output by the generator and the real clean image x MSE The calculation formula is as follows:

[0045]

[0046] Where h×w is the image space size and C is the number of channels; this loss term forces the generated image to be close to the real image in terms of pixel values;

[0047] S5.2: Use the structural similarity index to measure the structural similarity L between the generated image G(z) and the real image x SSIM , the formula is as follows:

[0048] L SSIM = 1-SSIM(G(z),x);

[0049] Among them, SSIM calculates the similarity of image brightness, contrast and structural information, with a value range of [0,1]. The calculation process includes image segmentation, Gaussian weighted mean estimation, variance and covariance calculation, and the formula is as follows:

[0050]

[0051] Among them, μ A 、μ B is the mean, σ A , σ B is the variance, σ AB is the covariance, C1, C2 are constants;

[0052] S5.4: Use the pre-trained visual feature extraction network to extract high-level semantic features of the generated image G(z) and the real image x, and calculate the mean square error in the feature space. Let φ(·) be the feature extraction function of the pre-trained network, and the output feature map dimension be h'×w'×c', where h', w' are the height and width of the feature map, and c' is the number of channels of the feature map. Then the perceptual loss L Perceptual The calculation formula is as follows:

[0053]

[0054] Among them, k∈[1,c′] represents the channel index of the feature map, which is used to traverse all channels of the feature map and calculate the eigenvalue difference of each channel at the spatial position (i,j);

[0055] S5.4: Weighted integration of MSE, SSIM, perceptual loss and Wasserstein loss to form the total loss function L of the generator G , the formula is as follows:

[0056] L G =λ MSE ·L MSE +λ SSIM ·L SSIM +λ Perceptual ·L Perceptual -E z～Pnoise [D(G(z))];

[0057] Among them, λ MSE ,λ SSIM ,λ Perceptual is a weighting parameter that balances the importance of different loss terms; E z～Pnoise [D(G(z))] is the inverse of the expected score of the discriminator on the generated sample, which corresponds to the goal of minimizing the Wasserstein distance of the generator in WGAN.

[0058] The step S6 further includes the following contents:

[0059] S6.1: Initialize the network parameters of the generator G and the discriminator D. Use a random normal distribution to initialize the weight matrix and the bias vector to 0 to break the symmetry and promote network learning. Set the training hyperparameters: the batch size is 32, which can be adjusted according to the computing resources, and the learning rate is 1×10 -4 , select the Adam optimizer and configure the parameters, and train for a number of iterations until the loss function converges; at the same time, clarify the components of the generator loss LG and the discriminator loss LD, including the MSE, SSIM, perceptual loss and Wasserstein distance terms defined in step S5, and determine the weighting parameters of each loss term to ensure synergy with the discriminator gradient penalty term in step S4;

[0060] S6.2: During each training round, perform the following operations:

[0061] Fixed generator, optimized discriminator: input batch real samples x and batch noisy samples z, the generator outputs denoised samples G(z); the discriminator scores x and G(z) and calculates the discriminator loss, the formula is as follows:

[0062] L D=-E[D(x)]+E[D(G(z))]+λ·GradientPenalty;

[0063] Among them, the first term in the formula represents the distribution of real images p real The discriminator score expectation of the sample in the image is , and the negative sign indicates that the discriminator goal is to maximize the true sample score; the second term represents the distribution of noisy images p noise The discriminator scores the expected value of the generated samples in the discriminator, and the discriminator reduces this value to distinguish the generated samples; the third term is the gradient penalty term, which is obtained by uniformly sampling the intermediate samples between the real samples and the generated samples. Calculate the discriminator output pair The expected square of the deviation of the L2 norm of the gradient from 1;

[0064] Fix the discriminator and optimize the generator: Keep the discriminator parameters unchanged, input the same batch of noisy samples z, and the generator output G(z); calculate the generator loss, the formula is shown in step S5.4;

[0065] The step S7 further includes the following contents:

[0066] S7.1: Preprocess the noisy low-dose CT image to be denoised in the same way as the training data in step S1. Then, based on the multi-energy spectrum characteristics of photon counting CT, separate the image into three channels: high-energy, low-energy, and full-energy, forming a three-dimensional input tensor [H noisy ,L noisy ,T noisy ], where H noisy represents a noisy high-energy image, L noisy represents a noisy low-energy image, T noisy represents a noisy full-scale image; this 3D input tensor is fed into the trained generator;

[0067] S7.2: The generator extracts features from the three-channel noisy images input separately; for the high-energy channel image H noisy , through 8 layers of convolutional kernel size of 5×5, the convolution operation is performed to extract the mid-scale features; for the low-energy channel image L noisy , using 8 convolution layers with a kernel size of 7×7 to extract the structural features of the global noise pattern; for the full channel image T noisy , use 8 layers of convolutional layers with a kernel size of 3×3 to extract texture and local detail features; after these operations, the deep feature maps F of each channel are obtained respectively Hnoisy 、F Lnoisy and F Tnoisy ;

[0068] S7.3: Extract the deep feature maps F of each channel Hnoisy 、F Lnoisy and F Tnoisy Input to the self-attention module; in the self-attention module, the feature map of each channel is globally averaged and pooled to generate a global description vector; then the similarity between these global description vectors is calculated, and the similarity results are normalized by the Softmax function to obtain the attention weights between each channel; finally, the feature maps of each channel are weighted fused according to these attention weights to obtain the multi-channel fusion feature map F fused ,as follows:

[0069]

[0070] S7.4: Multi-channel fusion feature map F fused Upsampling is performed using two deconvolution layers. The size of each deconvolution kernel is 4×4, the stride is 2, and the padding is 1. The size of the feature map is restored to the same size h×w as the input image. Finally, the Tanh activation function is used to output a denoised image with pixel values ​​in [0,1]. After upsampling, the feature map is converted into the final image, and a high-quality CT image with suppressed noise and preserved details is output.

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

[0072] 1. Efficient Utilization of Multi-Energy Spectral Information: This method extracts differentiated features from three independent input channels (a high-energy 5×5 core, a low-energy 7×7 core, and a full-energy 3×3 core). This, combined with a weighted fusion self-attention mechanism, addresses the inability of single-channel models to utilize information from multiple energy windows, significantly improving the ability to learn noise features. Experiments show that in walnut and mouse datasets, the peak signal-to-noise ratio (PSNR) of denoised images increases to 35.68 and 33.90, respectively. The structural similarity index (SSIM) improves by over 5% compared to traditional methods, significantly enhancing detail preservation.

[0073] 2. Optimized training stability and denoising: The WGAN architecture is adopted, and the discriminator is constrained by the Wasserstein distance and gradient penalty term, avoiding the mode collapse problem of traditional GANs and making the training process more stable. The generator loss function integrates the mean squared error (MSE), the structural similarity index (SSIM), and the perceptual loss to ensure that the denoised image is close to the real image in terms of numerical accuracy, structural similarity, and visual semantics, effectively solving the problems of edge blurring and texture loss.

[0074] 3. Multi-scale feature fusion and detail preservation: The multi-channel convolution kernel design (large kernels extract global noise, small kernels capture fine textures) combined with the self-attention mechanism enables the model to learn the generalized noise patterns of low-energy images while preserving the edge details of high-energy images and the texture information of all-energy images. This achieves the coordinated optimization of global and local noise patterns. The denoised image has clear details and the visual effect is significantly better than traditional methods such as BM3D and TV filtering. BRIEF DESCRIPTION OF THE DRAWINGS

[0075] Figure 1 Schematic diagram of the MCWGAN network generator structure of the photon counting CT image denoising method based on multi-channel WGAN of the present invention;

[0076] Figure 2 Schematic diagram of the MCWGAN network structure of the photon counting CT image denoising method based on multi-channel WGAN of the present invention;

[0077] Figure 3 This is a diagram showing the denoising effects of different methods of the photon counting CT image denoising method based on multi-channel WGAN in the present invention on the walnut dataset;

[0078] Figure 4 This is a comparison chart of the denoising effects of different methods of the photon counting CT image denoising method based on multi-channel WGAN in the mouse denoising task;

[0079] Figure 5 Schematic diagram of multi-channel feature fusion based on self-attention mechanism in the photon counting CT image denoising method based on multi-channel WGAN of the present invention;

[0080] Figure 6 Schematic diagram of the discriminator D of the photon counting CT image denoising method based on multi-channel WGAN in the present invention. DETAILED DESCRIPTION

[0081] The following will clearly and completely describe the technical solutions in the embodiments of the present invention in conjunction with the accompanying drawings. Obviously, the described embodiments are only part of the embodiments of the present invention, not all of the embodiments. Based on the embodiments of the present invention, all other embodiments obtained by ordinary technicians in this field without making creative efforts are within the scope of protection of the present invention.

[0082] Example: Figures 1-6 As shown, the present invention provides a technical solution.

[0083] The photon counting CT image denoising method based on multi-channel WGAN includes the following steps:

[0084] S1. Photon counting CT is used to collect high-energy, low-energy, and full-energy images. After noise reduction, correction, and normalization preprocessing, a three-channel input tensor is formed to provide a standardized data basis for subsequent feature extraction.

[0085] S2, extract deep feature maps from the three-channel images through independent convolution branches, and each layer is batch normalized and ReLU activated to obtain multi-scale semantic features;

[0086] S3. Perform global average pooling on the deep feature maps of each channel to generate a global vector, calculate the channel similarity weight through self-attention, and perform weighted fusion to form a fused feature map that focuses on the key features of multiple energy domains;

[0087] S4. Design a multi-layer convolutional neural network as the discriminator, use the Wasserstein distance to measure the distribution difference between the real and generated samples, and constrain it to meet the 1-Lipschitz condition through the gradient penalty term;

[0088] S5. Integrate Wasserstein loss, pixel error, structural similarity and perceptual loss to construct the total loss function of the generator to balance the numerical, structural and semantic consistency of the denoised image;

[0089] S6. Alternately train the generator and discriminator. Fix the generator to optimize the discriminator to improve the sample discrimination ability. Then fix the discriminator to optimize the generator to minimize the total loss. Iterate until the model converges.

[0090] S7. The preprocessed low-dose noisy CT image is input into the trained generator. After multi-channel feature extraction, self-attention fusion and upsampling, a high-quality image with noise suppression and detail preservation is output.

[0091] The step S1 further includes the following contents:

[0092] S1.1: The high-energy image H, low-energy image L, and full-energy image T of the same scanned object are acquired by a photon counting CT device, including the clean image NDCT obtained by normal-dose scanning and the noisy image LDCT obtained by low-dose scanning. Both are two-dimensional grayscale images with a uniform size of h×w, i.e., height×width. Among them, the clean image NDCT obtained by normal-dose scanning refers to the photon counting CT image obtained by standard radiation dose scanning, which is considered as a real image without significant noise after preprocessing; the noisy image LDCT obtained by low-dose scanning refers to the photon counting CT image obtained by low-dose radiation scanning. The image, including noise due to insufficient dose, serves as the noisy input data for the generator. NDCT stands for normal-dose photon-counting CT, and LDCT stands for low-dose CT. Based on the multi-energy spectral imaging characteristics of PCCT, high-energy images reflect the tissue's absorption characteristics of high-energy X-rays, low-energy images correspond to low-energy X-ray absorption characteristics, and full-energy images integrate information from the entire energy range. PCCT is a photon-counting CT device. The acquired images are then preprocessed, including removing background noise, correcting for detector response differences, and normalizing pixel values ​​to the neural network input range, which is [0,1].

[0093] S1.2: The preprocessed high-energy, low-energy, and full-energy images are used as a set of inputs to the generator, and combined with an optional random noise vector Z to form an input tensor with 3 channels; specifically expressed as [H, L, T, Z], where H∈R h×w , L∈R h×w , T∈R h×w , corresponding to the two-dimensional grayscale matrices of high-energy, low-energy, and full-energy images, respectively, h and w are the height and width of the image; the number of input channels of the generator is 3, which does not include channel Z, which is an optional additional input, and then directly receives the tensor as input; the generator will subsequently perform independent feature extraction on each channel, using a 5×5 convolution kernel for the high-energy channel, a 7×7 convolution kernel for the low-energy channel, and a 3×3 convolution kernel for the full-energy channel to ensure that the construction of the input tensor maintains the independence and correspondence of the three channels.

[0094] The step S2 further includes the following contents:

[0095] S2.1: In view of the different characteristic features of high-energy, low-energy and full-energy images, an independent convolutional neural network branch is designed for each channel, namely the high-energy image channel, low-energy image channel and full-energy image channel; for the high-energy image channel, a 5×5 convolution kernel is used, with a total of 8 convolution layers; the number of convolution kernels in the first layer is determined according to the number of input channels, initially set to 3 channels, and the number of convolution kernels in each subsequent layer is gradually increased according to the network design, with 64 convolution kernels in the first layer, and each subsequent layer increases by multiples, with the final number of kernels being 512, which is used to extract mid-scale features and balance noise and detail information; for the low-energy image channel, a 7×7 convolution kernel is used, also containing 8 convolution layers; used to extract generalized noise patterns and global structural features contained in low-energy images; for the full-energy image channel, a 3×3 convolution kernel is used, with 8 convolution layers; the specific operation process of each convolution layer is as follows: after the input image is calculated by the convolution kernel, the features are normalized by the batch normalization layer, and then nonlinear mapping is introduced by the ReLU activation function;

[0096] S2.2: The 8 convolutional layers of each channel extract features from the input image layer by layer. The front convolutional layer captures shallow features, and the back convolutional layer gradually integrates them to form deep semantic features. The 5×5 kernel of the high-energy channel focuses on extracting medium-scale noise-signal mixed features in the middle convolutional layer. The 7×7 kernel of the low-energy channel covers a large area in the shallow layer to capture the noise distribution law. The 3×3 kernel of the full-energy channel focuses on retaining fine textures in the deep layer. Finally, each channel outputs the corresponding deep feature map F. H 、F L 、F T Among them, F H 、F L 、F T They are represented as high-energy, low-energy and full-energy deep feature maps respectively.

[0097] The step S3 further includes the following contents:

[0098] S3.1: For each channel deep feature map F generated in step S2 H 、F L 、F T Perform global average pooling operations respectively; global average pooling operation is referred to as GAP operation later; these feature maps are all three-dimensional tensors with dimension R C×h×w , where C is the number of channels of the feature map, which is determined by the configuration of the generator convolutional layer, and h×w is the spatial size of the feature map, that is, the height × width of the image; the GAP operation compresses the three-dimensional tensor into a one-dimensional global description vector v∈R by performing channel-by-channel mean calculation on the spatial dimension of each feature map. C , the specific calculation is:

[0099] v H =GAP(FH )∈R C , v L =GAP(F L )∈R C , v T =GAP(F T )∈R C ;

[0100] in,

[0101]

[0102] c∈{1,2,...,C} is the channel index, that is, all spatial pixels of each channel are averaged to generate the global statistical features of the channel;

[0103] S3.2: Based on the one-dimensional global description vector v obtained in sub-step 3.1 H ,v L ,v T , calculate the attention weights between channels through the self-attention mechanism; transform the one-dimensional global description vector vH,vL,vT∈R C Through the fully connected layer mapping to the query matrix Qi and the key matrix Kj, i, j∈{H,L,T}, calculate the similarity matrix A∈R 3×3 , where A i,j =Q i ·K j , characterizes the pairwise correlation between the three channels; the query matrix Qi and the key matrix Kj, i, j∈{H,L,T} correspond to high energy, low energy, and full energy channels respectively; the dimension after mapping remains R C ; Then calculate the similarity between channel i and channel j by matrix dot product to obtain the similarity matrix A∈R 3×3 , where A i,j =Q i ·K j ; Then perform Softmax processing on the similarity matrix to generate normalized attention weights. The formula is as follows:

[0104]

[0105] Among them, α H , α L , α T are the attention weights of high-energy, low-energy, and full-energy channels, respectively, and α H +α L +α T =1, indicating the contribution of each channel in feature fusion;

[0106] S3.3: Use the normalized attention weight α calculated in substep 3.2 H , α L, α T , for each channel deep feature map F H 、F L 、F T Perform weighted fusion, and the dimensions of each feature map are R C×h×w , the fusion calculation formula is as follows:

[0107] F fused =α H ·F H +α L ·F L +α T ·F T ;

[0108] Among them, F fused ∈R C×h×w The pixel value of each channel of the fused feature map is calculated by the weighted sum of the corresponding positions of the original feature map; the contribution of each channel feature is dynamically adjusted through the attention weight, and the generator will focus on the key features of different energy domains: the generalized noise pattern extracted by the 7×7 convolution kernel in the low-energy channel, the edge details extracted by the 5×5 kernel in the high-energy channel, and the texture information extracted by the 3×3 kernel in the all-energy channel.

[0109] The step S4 further includes the following contents:

[0110] S4.1: Design the network architecture of the discriminator. Construct it as a 6-layer convolutional neural network. Each layer uses a 3×3 convolution kernel. The configuration of each layer is as follows:

[0111] For the convolution layer configuration, the first two layers are 64 convolution kernels with stride 1 and 64 convolution kernels with stride 2, respectively. The next two layers are 128 convolution kernels with stride 1 and 128 convolution kernels with stride 2, and the next two layers are 256 convolution kernels with stride 1 and 256 convolution kernels with stride 2. Each layer is followed by a LeakyReLU activation function. The convolution operation with stride 2 reduces the spatial size of the input image by half according to the preset range. At the same time, the number of channels is doubled layer by layer, and the low-level pixel features are gradually extracted to the high-level semantic features.

[0112] For the fully connected layer configuration, after 6 convolution layers, a 1024-dimensional fully connected layer and a 1-dimensional output layer are sequentially connected to output the score of the input sample. This structure does not use batch normalization and meets the Wasserstein GAN's 1-Lipschitz constraint on the discriminator.

[0113] S4.2: The discriminator receives a real clean image sample x and a denoised image sample G(z) output by the generator. Image sample x comes from a normal-dose photon-counting CT scan, a noise-free image acquired and preprocessed in step S1.1; denoised image sample G(z) comes from a noisy low-dose CT image input to the generator in step S1.2 and processed in steps S2-S3. The approximate Wasserstein distance between the two is then calculated. The first two terms of the loss function are as follows:

[0114] L(D)=-E x～Preal [D(x)]+E z～Pnoise [D(G(z))];

[0115] Among them, L(D) is the loss function of the discriminator, which is used to measure the approximate Wasserstein distance between the true distribution and the generated distribution. E is the mathematical expectation, which represents the average operation of the sample distribution. D(x) is the score of the discriminator on the true sample x. D(G(z)) is the score of the discriminator on the generated sample G(z). Preal is the distribution of the NDCT true image, and Pnoise is the distribution of the LDCT noisy image. By maximizing L(D), the discriminator is forced to improve the score of the true sample and reduce the score of the generated sample, thereby guiding the generator to adjust the parameters so that the distribution of G(z) approaches the true distribution p. real ;

[0116] Then constrain the discriminator to meet the 1-Lipschitz condition and uniformly sample intermediate samples along the line connecting the real sample and the generated sample Calculate the gradient penalty term as follows:

[0117]

[0118] Where λ is a constant weighting parameter; Represents the distribution of the middle sample Expected operation of ; Represents the discriminator output for the intermediate sample The gradient vector reflects the slope of the discriminator in the sample space; Represents the L2 norm of the gradient vector, which measures the length of the gradient; Represents the core of the penalty term, forcing the gradient norm to be close to 1 so that the discriminator satisfies the 1-Lipschitz constraint;

[0119] S4.3: Perform adversarial training, alternately optimizing the discriminator and generator. Specifically: input a batch of real samples x and generated samples G(z), calculate the total loss L(D), which includes the Wasserstein distance term and the gradient penalty term, and update the discriminator parameters through backpropagation. The goal is to maximize L(D) and improve its ability to distinguish between real and generated samples; then rely on the gradient penalty term to constrain the discriminator gradient norm.

[0120] The step S5 further includes the following contents:

[0121] S5.1: Pixel-wise error L between the denoised image G(z) output by the generator and the real clean image x MSE The calculation formula is as follows:

[0122]

[0123] Where h×w is the image space size and C is the number of channels; this loss term forces the generated image to be close to the real image in terms of pixel values;

[0124] S5.2: Use the structural similarity index to measure the structural similarity L between the generated image G(z) and the real image x SSIM , the formula is as follows:

[0125] L SSIM = 1-SSIM(G(z),x);

[0126] Among them, SSIM calculates the similarity of image brightness, contrast and structural information, with a value range of [0,1]. The calculation process includes image segmentation, Gaussian weighted mean estimation, variance and covariance calculation, and the formula is as follows:

[0127]

[0128] Among them, μ A 、μ B is the mean, σ A , σ B is the variance, σ AB is the covariance, C1, C2 are constants;

[0129] S5.4: Use the pre-trained visual feature extraction network to extract high-level semantic features of the generated image G(z) and the real image x, and calculate the mean square error in the feature space. Let φ(·) be the feature extraction function of the pre-trained network, and the output feature map dimension be h'×w'×c', where h', w' are the height and width of the feature map, and c' is the number of channels of the feature map. Then the perceptual loss L Perceptual The calculation formula is as follows:

[0130]

[0131] Among them, k∈[1,c′] represents the channel index of the feature map, which is used to traverse all channels of the feature map and calculate the eigenvalue difference of each channel at the spatial position (i,j);

[0132] S5.4: Weighted integration of MSE, SSIM, perceptual loss and Wasserstein loss to form the total loss function L of the generatorG , the formula is as follows:

[0133] L G =λ MSE·LMSE +λ SSIM ·L SSIM +λ Perceptual ·L Perceptual -E z～Pnoise [D(G(z))];

[0134] Among them, λ MSE ,λ SSIM ,λ Perceptual is a weighting parameter that balances the importance of different loss terms; E z～Pnoise [D(G(z))] is the inverse of the expected score of the discriminator on the generated sample, which corresponds to the goal of minimizing the Wasserstein distance of the generator in WGAN.

[0135] The step S6 further includes the following contents:

[0136] S6.1: Initialize the network parameters of the generator G and the discriminator D. Use a random normal distribution to initialize the weight matrix and the bias vector to 0 to break the symmetry and promote network learning. Set the training hyperparameters: the batch size is 32, which can be adjusted according to the computing resources, and the learning rate is 1×10 -4 , select the Adam optimizer and configure the parameters, and train for a number of iterations until the loss function converges; at the same time, clarify the components of the generator loss LG and the discriminator loss LD, including the MSE, SSIM, perceptual loss and Wasserstein distance terms defined in step S5, and determine the weighting parameters of each loss term to ensure synergy with the discriminator gradient penalty term in step S4;

[0137] S6.2: During each training round, perform the following operations:

[0138] Fixed generator, optimized discriminator: input batch real samples x and batch noisy samples z, the generator outputs denoised samples G(z); the discriminator scores x and G(z) and calculates the discriminator loss, the formula is as follows:

[0139] L D =-E[D(x)]+E[D(G(z))]+λ·GradientPenalty;

[0140] Among them, the first term in the formula represents the distribution of real images p real The discriminator score expectation of the sample in the image is , and the negative sign indicates that the discriminator goal is to maximize the true sample score; the second term represents the distribution of noisy images p noiseThe discriminator scores the expected value of the generated samples in the discriminator, and the discriminator reduces this value to distinguish the generated samples; the third term is the gradient penalty term, which is obtained by uniformly sampling the intermediate samples between the real samples and the generated samples. Calculate the discriminator output pair The expected square of the deviation of the L2 norm of the gradient from 1;

[0141] Fix the discriminator and optimize the generator: Keep the discriminator parameters unchanged, input the same batch of noisy samples z, and the generator output G(z); calculate the generator loss, the formula is shown in step S5.4;

[0142] The step S7 further includes the following contents:

[0143] S7.1: Preprocess the noisy low-dose CT image to be denoised in the same way as the training data in step S1. Then, based on the multi-energy spectrum characteristics of photon counting CT, separate the image into three channels: high-energy, low-energy, and full-energy, forming a three-dimensional input tensor [H noisy ,L noisy ,T noisy ], where H noisy represents a noisy high-energy image, L noisy represents a noisy low-energy image, T noisy represents a noisy full-scale image; this 3D input tensor is fed into the trained generator;

[0144] S7.2: The generator extracts features from the three-channel noisy images input separately; for the high-energy channel image H noisy , through 8 layers of convolutional kernel size of 5×5, the convolution operation is performed to extract the mid-scale features; for the low-energy channel image L noisy , using 8 convolution layers with a kernel size of 7×7 to extract the structural features of the global noise pattern; for the full channel image T noisy , use 8 layers of convolutional layers with a kernel size of 3×3 to extract texture and local detail features; after these operations, the deep feature maps F of each channel are obtained respectively Hnoisy 、F Lnoisy and F Tnoisy ;

[0145] S7.3: Extract the deep feature maps F of each channel Hnoisy 、F Lnoisy and F TnoisyInput to the self-attention module; in the self-attention module, the feature map of each channel is globally averaged and pooled to generate a global description vector; then the similarity between these global description vectors is calculated, and the similarity results are normalized by the Softmax function to obtain the attention weights between each channel; finally, the feature maps of each channel are weighted fused according to these attention weights to obtain the multi-channel fusion feature map F fused ,as follows:

[0146]

[0147] S7.4: Multi-channel fusion feature map F fused Upsampling is performed using two deconvolution layers. The size of each deconvolution kernel is 4×4, the stride is 2, and the padding is 1. The size of the feature map is restored to the same size h×w as the input image. Finally, the Tanh activation function is used to output a denoised image with pixel values ​​in [0,1]. After upsampling, the feature map is converted into the final image, and a high-quality CT image with suppressed noise and preserved details is output.

[0148] This example uses a Siemens NAEOTOM Alpha photon-counting CT system (equipped with a cadmium telluride detector, an energy resolution of 1 keV, and a spatial resolution of 0.25 mm) to perform whole-body scans on patients with multiple myeloma. Three energy spectrum acquisition modes are set: high energy 140 keV, low energy 100 keV, and full energy mode. The radiation dose is controlled at 1.61 mSv (CTDIvol), which is only 13.5% of the traditional CT (11.9 mSv). After acquisition, the original image undergoes three steps of preprocessing:

[0149] Noise suppression: A 3×3 Gaussian filter (σ=1.5) was used to remove high-frequency noise, combined with a bilateral filter (d=30, sigmaColor=35, sigmaSpace=15) to preserve edge details.

[0150] Response correction: Corrects detector response differences based on air calibration tables, eliminates hardware non-uniformity, and improves detector signal consistency to 99.2%;

[0151] Normalization: Map the image grayscale values ​​to the range [0, 1] to form an input tensor of size 256×256×3, which serves as the standardized input for subsequent models.

[0152] The generator uses a three-branch convolutional neural network, each consisting of four convolutional layers (kernel sizes of 5×5, 3×3, 3×3, and 3×3, with strides of 2, 1, 1, and 1, respectively). Each layer is followed by batch normalization and Reluctant Unit (ReLU) activation to extract multi-scale features from different energy domains. A self-attention module performs global average pooling on the deep feature maps of each branch, generating a 1×1×512 global description vector. Weights are generated by calculating the inter-channel similarity matrix, and weighted fusion is performed to form a fused feature map that focuses on key features. The upsampling module restores the feature map to a 256×256 resolution through deconvolution (kernel size 4×4, stride 2) and bilinear interpolation, outputting a denoised CT image. The discriminator is a four-layer convolutional neural network (kernel size 3×3, stride 2). Each layer is followed by a LeakyReLU (α=0.2) activation function, ultimately outputting a single scalar value representing the probability of the sample being genuine. The Wasserstein distance is used to measure the distribution difference between real and generated samples, and a gradient penalty term (weight 10) is introduced to constrain the discriminator to meet the 1-Lipschitz condition. The loss function combines the Wasserstein loss (weight 1), mean squared error (MSE, weight 0.5), structural similarity (SSIM, weight 0.3), and VGG16 perceptual loss (weight 0.2). The generator and discriminator are trained alternately using the Adam optimizer (learning rate 2e-4, β1 = 0.5), with a batch size of 64 and 1000 iterations until convergence.

[0153] The noisy low-dose CT image (noise level σ = 30) is input into the trained generator, and the output image is evaluated:

[0154] Objective indicators: Peak signal-to-noise ratio (PSNR) increased by 3.2dB, structural similarity index (SSIM) increased from 0.82 to 0.95, and iodine signal-to-noise ratio increased by 40%;

[0155] Subjective evaluation: Professionals rated the generated images on a 5-point scale (1 = poor, 5 = excellent) significantly higher than the traditional filtering method in terms of detection of osteolytic lesions, vascular clarity, and soft tissue contrast (4.2 ± 0.3 vs 3.1 ± 0.5).

[0156] It will be apparent to those skilled in the art that the present invention is not limited to the details of the exemplary embodiments described above and that the invention can be embodied in other specific forms without departing from the spirit or essential characteristics of the invention. Therefore, the embodiments should be considered in all respects as illustrative and non-restrictive, and the scope of the invention is defined by the appended claims, not the foregoing description, and all variations within the meaning and range of equivalents of the claims are intended to be included therein. Any reference sign in a claim should not be construed as limiting the claim to which it relates.

Claims

1. A photon counting CT image denoising method based on a multi-channel WGAN, characterized by: The following steps are involved: S1. Photon counting CT is used to collect high-energy, low-energy, and full-energy images. After noise reduction, correction, and normalization preprocessing, a three-channel input tensor is formed to provide a standardized data basis for subsequent feature extraction. S2, extract deep feature maps from the three-channel images through independent convolution branches, and each layer is batch normalized and ReLU activated to obtain multi-scale semantic features; S3. Perform global average pooling on the deep feature maps of each channel to generate a global vector, calculate the channel similarity weight through self-attention, and perform weighted fusion to form a fused feature map that focuses on the key features of multiple energy domains; S4. Design a multi-layer convolutional neural network as the discriminator, use the Wasserstein distance to measure the distribution difference between the real and generated samples, and constrain it to meet the 1-Lipschitz condition through the gradient penalty term; S5. Integrate Wasserstein loss, pixel error, structural similarity and perceptual loss to construct the total loss function of the generator to balance the numerical, structural and semantic consistency of the denoised image; S6. Alternately train the generator and discriminator. Fix the generator to optimize the discriminator to improve the sample discrimination ability. Then fix the discriminator to optimize the generator to minimize the total loss. Iterate until the model converges. S7. The preprocessed low-dose noisy CT image is input into the trained generator. After multi-channel feature extraction, self-attention fusion and upsampling, a high-quality image with noise suppression and detail preservation is output.

2. The photon counting CT image denoising method based on multi-channel WGAN according to claim 1 is characterized in that: The step S1 further includes the following contents: S1.1: The high-energy image H, low-energy image L, and full-energy image T of the same scanned object are acquired by a photon counting CT device, including the clean image NDCT obtained by normal-dose scanning and the noisy image LDCT obtained by low-dose scanning. Both are two-dimensional grayscale images with a uniform size of h×w, i.e., height×width. Among them, the clean image NDCT obtained by normal-dose scanning refers to the photon counting CT image obtained by standard radiation dose scanning, which is considered as a real image without significant noise after preprocessing; the noisy image LDCT obtained by low-dose scanning refers to the photon counting CT image obtained by low-dose radiation scanning. The image, including noise due to insufficient dose, serves as the noisy input data for the generator. NDCT stands for normal-dose photon-counting CT, and LDCT stands for low-dose CT. Based on the multi-energy spectral imaging characteristics of PCCT, high-energy images reflect the tissue's absorption characteristics of high-energy X-rays, low-energy images correspond to low-energy X-ray absorption characteristics, and full-energy images integrate information from the entire energy range. PCCT is a photon-counting CT device. The acquired images are then preprocessed, including removing background noise, correcting for detector response differences, and normalizing pixel values ​​to the neural network input range, which is [0,1]. S1.2: The preprocessed high-energy, low-energy, and full-energy images are used as a set of inputs to the generator, and combined with an optional random noise vector Z to form an input tensor with 3 channels; specifically expressed as [H, L, T, Z], where H∈R h×w , L∈R h×w , T∈R h ×w , corresponding to the two-dimensional grayscale matrices of high-energy, low-energy, and full-energy images, respectively, h and w are the height and width of the image; the number of input channels of the generator is 3, which does not include channel Z, which is an optional additional input, and then directly receives the tensor as input; the generator will subsequently perform independent feature extraction on each channel, using a 5×5 convolution kernel for the high-energy channel, a 7×7 convolution kernel for the low-energy channel, and a 3×3 convolution kernel for the full-energy channel to ensure that the construction of the input tensor maintains the independence and correspondence of the three channels.

3. The photon counting CT image denoising method based on multi-channel WGAN according to claim 1, characterized in that: The step S2 further includes the following contents: S2.1: In view of the different characteristic features of high-energy, low-energy and full-energy images, an independent convolutional neural network branch is designed for each channel, namely the high-energy image channel, low-energy image channel and full-energy image channel; for the high-energy image channel, a 5×5 convolution kernel is used, with a total of 8 convolution layers; the number of convolution kernels in the first layer is determined according to the number of input channels, initially set to 3 channels, and the number of convolution kernels in each subsequent layer is gradually increased according to the network design, with 64 convolution kernels in the first layer, and each subsequent layer increases by multiples, with the final number of kernels being 512, which is used to extract mid-scale features and balance noise and detail information; for the low-energy image channel, a 7×7 convolution kernel is used, also containing 8 convolution layers; used to extract generalized noise patterns and global structural features contained in low-energy images; for the full-energy image channel, a 3×3 convolution kernel is used, with 8 convolution layers; the specific operation process of each convolution layer is as follows: after the input image is calculated by the convolution kernel, the features are normalized by the batch normalization layer, and then nonlinear mapping is introduced by the ReLU activation function; S2.2: The 8 convolutional layers of each channel extract features from the input image layer by layer. The front convolutional layer captures shallow features, and the back convolutional layer gradually integrates them to form deep semantic features. The 5×5 kernel of the high-energy channel focuses on extracting medium-scale noise-signal mixed features in the middle convolutional layer. The 7×7 kernel of the low-energy channel covers a large area in the shallow layer to capture the noise distribution law. The 3×3 kernel of the full-energy channel focuses on retaining fine textures in the deep layer. Finally, each channel outputs the corresponding deep feature map F. H 、F L 、F T Among them, F H 、F L 、F T They are represented as high-energy, low-energy and full-energy deep feature maps respectively.

4. The photon counting CT image denoising method based on multi-channel WGAN according to claim 3 is characterized by: The step S3 further includes the following contents: S3.1: For each channel deep feature map F generated in step S2 H 、F L 、F T Perform global average pooling operations respectively; global average pooling operation is referred to as GAP operation later; these feature maps are all three-dimensional tensors with dimension R C×h×w , where C is the number of channels of the feature map, which is determined by the configuration of the generator convolutional layer, and h×w is the spatial size of the feature map, that is, the height × width of the image; the GAP operation compresses the three-dimensional tensor into a one-dimensional global description vector v∈R by performing channel-by-channel mean calculation on the spatial dimension of each feature map. C , the specific calculation is: v H =GAP(F H )∈R C ,v L =GAP(F L )∈R C ,v T =GAP(F T )∈R C ; in, c∈{1,2,...,C} is the channel index, that is, all spatial pixels of each channel are averaged to generate the global statistical features of the channel; S3.2: Based on the one-dimensional global description vector v obtained in sub-step 3.1 H ,v L ,v T , calculate the attention weights between channels through the self-attention mechanism; transform the one-dimensional global description vector vH,vL,vT∈R C Through the fully connected layer mapping to the query matrix Qi and the key matrix Kj, i, j∈{H,L,T}, calculate the similarity matrix A∈R 3×3 , where A i,j =Q i ·K j , characterizes the pairwise correlation between the three channels; the query matrix Qi and the key matrix Kj, i, j∈{H,L,T} correspond to high energy, low energy, and full energy channels respectively; the dimension after mapping remains R C ; Then calculate the similarity between channel i and channel j by matrix dot product to obtain the similarity matrix A∈R 3×3 , where A i,j =Q i ·K j ; Then perform Softmax processing on the similarity matrix to generate normalized attention weights. The formula is as follows: Among them, α H , α L , α T are the attention weights of high-energy, low-energy, and full-energy channels, respectively, and α H +α L +α T =1, indicating the contribution of each channel in feature fusion; S3.3: Use the normalized attention weight α calculated in substep 3.2 H , α L , α T , for each channel deep feature map F H 、F L 、F T Perform weighted fusion, and the dimensions of each feature map are R C×h×w , the fusion calculation formula is as follows: F fused =a H ·F H +a L ·F L +a T ·F T ; Among them, F fused ∈R C×h×w The pixel value of each channel of the fused feature map is calculated by the weighted sum of the corresponding positions of the original feature map; the contribution of each channel feature is dynamically adjusted through the attention weight, and the generator will focus on the key features of different energy domains: the generalized noise pattern extracted by the 7×7 convolution kernel in the low-energy channel, the edge details extracted by the 5×5 kernel in the high-energy channel, and the texture information extracted by the 3×3 kernel in the all-energy channel.

5. The photon counting CT image denoising method based on multi-channel WGAN according to claim 2, characterized in that: The step S4 further includes the following contents: S4.1: Design the network architecture of the discriminator. Construct it as a 6-layer convolutional neural network. Each layer uses a 3×3 convolution kernel. The configuration of each layer is as follows: For the convolution layer configuration, the first two layers are 64 convolution kernels with stride 1 and 64 convolution kernels with stride 2, respectively. The next two layers are 128 convolution kernels with stride 1 and 128 convolution kernels with stride 2, and the next two layers are 256 convolution kernels with stride 1 and 256 convolution kernels with stride 2. Each layer is followed by a LeakyReLU activation function. The convolution operation with stride 2 reduces the spatial size of the input image by half according to the preset range. At the same time, the number of channels is doubled layer by layer, and the low-level pixel features are gradually extracted to the high-level semantic features. For the fully connected layer configuration, after 6 convolution layers, a 1024-dimensional fully connected layer and a 1-dimensional output layer are sequentially connected to output the score of the input sample; This structure does not use batch normalization and meets the 1-Lipschitz constraint of Wasserstein GAN on the discriminator; S4.2: The discriminator receives a real clean image sample x and a denoised image sample G(z) output by the generator. Image sample x comes from a normal-dose photon-counting CT scan, a noise-free image acquired and preprocessed in step S1.1; denoised image sample G(z) comes from a noisy low-dose CT image input to the generator in step S1.2 and processed in steps S2-S3. The approximate Wasserstein distance between the two is then calculated. The first two terms of the loss function are as follows: L(D)=-E x～Preal [D(x)]+E z ~P noise [D(G(z))]; Among them, L(D) is the loss function of the discriminator, which is used to measure the approximate Wasserstein distance between the true distribution and the generated distribution. E is the mathematical expectation, which represents the average operation of the sample distribution. D(x) is the score of the discriminator on the true sample x. D(G(z)) is the score of the discriminator on the generated sample G(z). Preal is the distribution of the NDCT true image, and Pnoise is the distribution of the LDCT noisy image. By maximizing L(D), the discriminator is forced to improve the score of the true sample and reduce the score of the generated sample, thereby guiding the generator to adjust the parameters so that the distribution of G(z) approaches the true distribution p. real ; Then constrain the discriminator to meet the 1-Lipschitz condition and uniformly sample intermediate samples along the line connecting the real sample and the generated sample Calculate the gradient penalty term as follows: Where λ is a constant weighting parameter; Represents the distribution of the middle sample Expected operation of ; Represents the discriminator output for the intermediate sample The gradient vector reflects the slope of the discriminator in the sample space; Represents the L2 norm of the gradient vector, which measures the length of the gradient; Represents the core of the penalty term, forcing the gradient norm to be close to 1 so that the discriminator satisfies the 1-Lipschitz constraint; S4.3: Perform adversarial training, alternately optimizing the discriminator and generator. Specifically: input a batch of real samples x and generated samples G(z), calculate the total loss L(D), which includes the Wasserstein distance term and the gradient penalty term, and update the discriminator parameters through backpropagation. The goal is to maximize L(D) and improve its ability to distinguish between real and generated samples; then rely on the gradient penalty term to constrain the discriminator gradient norm.

6. The photon counting CT image denoising method based on multi-channel WGAN according to claim 1, characterized in that: The step S5 further includes the following contents: S5.1: Pixel-wise error L between the denoised image G(z) output by the generator and the real clean image x MSE The calculation formula is as follows: Where h×w is the image space size and C is the number of channels; this loss term forces the generated image to be close to the real image in terms of pixel values; S5.2: Use the structural similarity index to measure the structural similarity L between the generated image G(z) and the real image x SSIM , the formula is as follows: L SSIM =1-SSIM(G(z),x); Among them, SSIM calculates the similarity of image brightness, contrast and structural information, with a value range of [0,1]. The calculation process includes image segmentation, Gaussian weighted mean estimation, variance and covariance calculation, and the formula is as follows: Among them, μ A 、μ B is the mean, σ A , σ B is the variance, σ AB is the covariance, C1, C2 are constants; S5.4: Use the pre-trained visual feature extraction network to extract high-level semantic features of the generated image G(z) and the real image x, and calculate the mean square error in the feature space. Let φ(·) be the feature extraction function of the pre-trained network, and the output feature map dimension be h'×w'×c', where h', w' are the height and width of the feature map, and c' is the number of channels of the feature map. Then the perceptual loss L Perceptual The calculation formula is as follows: Among them, k∈[1,c′] represents the channel index of the feature map, which is used to traverse all channels of the feature map and calculate the eigenvalue difference of each channel at the spatial position (i,j); S5.4: Weighted integration of MSE, SSIM, perceptual loss and Wasserstein loss to form the total loss function L of the generator G , the formula is as follows: I G =λMSE·L MSE +λ SSIM ·L SSIM +λ Perceptual ·L Perceptual - E z~P noise [D(G(z))]; Among them, λ MSE ,λ SSIM ,λ Perceptual is a weighting parameter that balances the importance of different loss terms; E z～Pnoise [D(G(z))] is the inverse of the expected score of the discriminator on the generated sample, which corresponds to the goal of minimizing the Wasserstein distance of the generator in WGAN.

7. The photon counting CT image denoising method based on multi-channel WGAN according to claim 6, characterized in that: The step S6 further includes the following contents: S6.1: Initialize the network parameters of the generator G and the discriminator D. Use a random normal distribution to initialize the weight matrix and the bias vector to 0 to break the symmetry and promote network learning. Set training hyperparameters: batch size is 32, which can be adjusted according to computing resources, and learning rate is 1×10 -4 , select the Adam optimizer and configure the parameters, and train for a number of iterations until the loss function converges; at the same time, clarify the components of the generator loss LG and the discriminator loss LD, including the MSE, SSIM, perceptual loss and Wasserstein distance terms defined in step S5, and determine the weighting parameters of each loss term to ensure synergy with the discriminator gradient penalty term in step S4; S6.2: During each training round, perform the following operations: Fixed generator, optimized discriminator: input batch real samples x and batch noisy samples z, the generator outputs denoised samples G(z); the discriminator scores x and G(z) and calculates the discriminator loss, the formula is as follows: L D =-E[D(x)]+E[D(G(z))]+λ·GradientPenalty; Among them, the first term in the formula represents the distribution of real images p real The discriminator score expectation of the sample in the image is , and the negative sign indicates that the discriminator goal is to maximize the true sample score; the second term represents the distribution of noisy images p noise The discriminator scores the expected value of the generated samples in the discriminator, and the discriminator reduces this value to distinguish the generated samples; the third term is the gradient penalty term, which is obtained by uniformly sampling the intermediate samples between the real samples and the generated samples. Calculate the discriminator output pair The expected square of the deviation of the L2 norm of the gradient from 1; Fix the discriminator and optimize the generator: Keep the discriminator parameters unchanged, input the same batch of noisy samples z, and the generator output G(z); calculate the generator loss, the formula is shown in step S5.

4.

8. The photon counting CT image denoising method based on multi-channel WGAN according to claim 1, characterized in that: The step S7 further includes the following contents: S7.1: Preprocess the noisy low-dose CT image to be denoised. This preprocessing method is consistent with the preprocessing method of the training data in step S1. Then, according to the multi-energy spectrum characteristics of photon counting CT, the image is separated into three channels: high energy, low energy and full energy, forming a three-dimensional input tensor [H noisy ,L noisy ,T noisy ], where H noisy represents a noisy high-energy image, L noisy represents a noisy low-energy image, T noisy represents a noisy full-scale image; this 3D input tensor is fed into the trained generator; S7.2: The generator extracts features from the three-channel noisy images input separately; for the high-energy channel image H noisy , convolution operation is performed through 8 layers of convolutional kernels with a size of 5×5 to extract the mid-scale features; For the low-energy channel image L noisy , using 8 convolution layers with a kernel size of 7×7 to extract the structural features of the global noise pattern; for the full channel image T noisy , use 8 layers of convolutional layers with a kernel size of 3×3 to extract texture and local detail features; after these operations, the deep feature maps F of each channel are obtained respectively Hnoisy 、F Lnoisy and F Tnoisy ; S7.3: Extract the deep feature maps F of each channel Hnoisy 、F Lnoisy and F Tnoisy Input to the self-attention module; in the self-attention module, the feature map of each channel is globally averaged and pooled to generate a global description vector; then the similarity between these global description vectors is calculated, and the similarity results are normalized by the Softmax function to obtain the attention weights between each channel; finally, the feature maps of each channel are weighted fused according to these attention weights to obtain the multi-channel fusion feature map F fused ,as follows: S7.4: Multi-channel fusion feature map F fused Upsampling is performed using two deconvolution layers. The size of each deconvolution kernel is 4×4, the stride is 2, and the padding is 1. The size of the feature map is restored to the same size h×w as the input image. Finally, the Tanh activation function is used to output a denoised image with pixel values ​​in [0,1]. After upsampling, the feature map is converted into the final image, and a high-quality CT image with suppressed noise and preserved details is output.

Images