Unsupervised hyperspectral image super-resolution method based on matrix factorization network

Abstract

The application discloses an unsupervised hyperspectral image super-resolution method based on a matrix decomposition network and belongs to the technical field of hyperspectral image processing. The application is used for processing a hyperspectral image, generating a simulated low spatial resolution hyperspectral image Y and a high spatial resolution multispectral image Z; first, inputting the generated data pair (Y, Z) into a designed auto-encoder network, training iteration to obtain a point spread function and a spectral response function; for a target high spatial resolution hyperspectral image X, the target high spatial resolution hyperspectral image X can be assumed to be a linear combination of a terminal member matrix A and a corresponding abundance matrix S, that is, X = AS, a spectral and spatial degradation model is combined to model, a deep CP decomposition module is designed to calculate A, A and S are iteratively solved, and finally a fusion result is obtained. The application can obtain more rich spectral and spatial features, obtain a better fusion result, and has good performance in practice.

Description

Technical Field

[0001] This invention relates to an unsupervised hyperspectral image super-resolution method based on matrix factorization networks, belonging to the field of hyperspectral image processing technology. Background Technology

[0002] Hyperspectral remote sensing, also known as high-spectral-resolution remote sensing, utilizes imaging techniques and spectral features to acquire a series of continuous, narrow-band image data with high spectral resolution from objects of interest. Hyperspectral images not only provide spatial information about land cover targets but also offer rich spectral information. They have been widely applied in various fields such as geological exploration, climate monitoring, crop detection, and national defense.

[0003] Hyperspectral image classification is an important task in the field of remote sensing applications. In recent years, with the continuous development of deep learning, hyperspectral classification methods based on deep networks have also been proven to be one of the effective means of classification. Supervised deep learning methods can be regarded as a nonlinear mapping from the feature space to the label space, which can extract higher-level abstract semantic features from the original data.

[0004] Due to the limited incident energy, there is always a trade-off between spectral resolution, spatial resolution, and signal-to-noise ratio when designing imaging sensors. Therefore, spatial resolution of HSI is often sacrificed, hindering subsequent tasks. Conversely, traditional multispectral imaging (MSI) can achieve higher spatial resolution at lower spectral resolution. An economical HSI super-resolution solution is to record low spatial resolution HSI (LR-HSI) and high spatial resolution MSI (HR-MSI) and fuse them into a target high spatial resolution HSI (HR-HSI). Hyperspectral imaging super-resolution solutions overcome the limitations of imaging sensors and have demonstrated good performance in practice.

[0005] Hyperspectral remote sensing images (HSI) typically have high spectral resolution and low spatial resolution. Conversely, multispectral images (MSI) typically have low spectral resolution and high spatial resolution. The image inference problem that combines the hyperspectral and high spatial resolution of HSI and MSI is a data fusion problem, namely hyperspectral image super-resolution. Hyperspectral images contain rich spectral details, which are beneficial for various computer vision tasks such as face recognition, object tracking, and modeling.

[0006] Existing hyperspectral image super-resolution methods mainly include component replacement (CS), multi-resolution analysis (MAR), model-driven methods, and deep learning-based methods. However, these methods cannot simultaneously acquire rich spectral and spatial features, resulting in poor final fusion results. Summary of the Invention

[0007] To address the challenge of obtaining more spatial and spectral features and improving fusion results, this invention provides a hyperspectral image super-resolution method, comprising:

[0008] Step 1: Process the hyperspectral image to generate a simulated low spatial resolution hyperspectral image Y and a high spatial resolution multispectral image Z as training data;

[0009] Step 2: Input the data pair (Y, Z) generated in Step 1 into the autoencoder network, and train iteratively to obtain the point spread function B and the spectral response function R;

[0010] Step 3: For the target image, the high spatial resolution hyperspectral image X, it is assumed to be a linear combination of the terminal member matrix and the corresponding abundance matrix; the initial terminal member matrix A0 and the initial abundance matrix S0 are used as trainable parameters, and the estimated B and R in Step 2 are input into the deep CP decomposition network to solve for the desired terminal member matrix A through the deep CP decomposition network.

[0011] Step 4: Combine the point spread function B and spectral response function R obtained in Step 2, and combine the degradation model and matrix factorization model to model the target, design regularization terms, and iteratively solve the terminal member matrix A and the corresponding abundance matrix S through network expansion to finally obtain the fused image.

[0012] Optionally, step 2 specifically includes:

[0013] Step 2.1: The degradation models for the low spatial resolution hyperspectral image Y and the high spatial resolution multispectral image Z are respectively modeled as follows:

[0014] Z = RX + N m Y = XB + N h

[0015] Where, N m and N h These represent different signal-to-noise ratios, and X represents a high spatial resolution hyperspectral image of the target.

[0016] Combining the two equations above, we get: ZB = RY, and the final objective function is: min||ZB-RY||;

[0017] Step 2.2: Solve using an autoencoder network:

[0018] Input the generated data pair (Y, Z), randomly initialize the point spread function B and the spectral response function R. The width and height of B are both (2*ratio)-1, where ratio is the sampling rate, i.e. the ratio of the spatial resolution of the two data points. The size of R is (C, c), where C represents the number of bands in the hyperspectral data and c represents the number of bands in the multispectral data. Use B and R as trainable parameters of the network.

[0019] Step 2.3: Use R as the convolution kernel to implement tensor multiplication in the form of convolution. Similarly, apply it to ZB to optimize the objective function min||ZB-RY||. Iterate to solve Z1 and Y1, where Z1 and Y1 are obtained by the second-order degradation of the learned B and R. Normalize the obtained Z1 and Y1 and save B and R.

[0020] Optionally, step 3 specifically includes:

[0021] Step 3.1: Based on the matrix factorization model: X = AS, establish the objective function: min||X-AS||, where A represents the terminal member matrix and S represents the corresponding abundance matrix;

[0022] Step 3.2: Randomly initialize the terminal member matrix A and the abundance matrix S to obtain the initial terminal member matrix A0 and the corresponding initial abundance matrix S0, and obtain the initial target image X0 = A0S0; where A has a size of (1, rank, H, W), where rank is the set CP rank, H and W are the length and width of the multispectral data, and S has a size of (C, rank, 1, 1), where C is the number of bands of the hyperspectral data. The point spread function B and the spectral response function R obtained in Step 2 are used as inputs.

[0023] Step 3.3: Input the initial terminal member matrix A0 into the low-rank tensor generation module LRTG. First, extract shallow features through a convolutional Conv, and then use Conv-ReLU-Conv to generate the input feature map. The size of the spatial kernel is 3x3. Input the extracted feature map into the depth CP decomposition module. First, perform global pooling GAP on the channel, height and width dimensions respectively to obtain vectors C1, C2 and C3. Perform Kroneckor product on the three to generate tensor O. Then, introduce residual learning to learn tensors Oi (i = 1, 2... r) of different frequencies:

[0024] For the generated O i Extracting the residual portion, treating it as O i Unrecoverable high-frequency information; the remaining portion is processed again by the low-rank tensor generation module LRTG to generate O. i+1(i = 1, 2, ..., r); Finally, the low-rank tensor learning module LRTL is entered. Given a predefined CP rank r, r tensors {O1, O2, ..., O} are generated through r iterations of the low-rank tensor generation module. r Finally, all generated tensors are aggregated into a low-rank tensor, and convolutional layers are used to aggregate the stacked tensors into the desired end member matrix A.

[0025] Optionally, step 4 includes:

[0026] From the degeneracy model Z = RX + N m Y = XB + N h The optimization objective is obtained by modeling the simultaneous matrix decomposition formula X = AS: min||Z-RAS||1+||Y-ASB||1, and the design regularization term ||A-A0||1 is added. Using the point spread function B and spectral response function R estimated in step 2, as well as the initialized A0 and S0, A and S are solved:

[0027] Solve for A i Solve using a deep CP decomposition network, with input A. i-1 Output A i i = 1, 2, ..., n, where n is the iteration number, which is the number of subsequent network training iterations;

[0028] Solve for S i S is set as a learnable parameter that is automatically learned and updated as the network is trained.

[0029] Iterate through the above steps until convergence, solve for A and S, and finally obtain the fused result X = AS.

[0030] Optionally, the loss function for training the deep CP decomposition network is as follows:

[0031] Loss=||Z2-Z||1+||Y2-Y||1+a||A-A0||1

[0032] Z2 and Y2 are low-spatial-resolution hyperspectral and high-spatial-resolution multispectral images obtained by further degrading the generated X. The degradation uses real B and R, and parameter a is a constant.

[0033] Optionally, parameter a can be set to 0.008.

[0034] Optionally, it also includes: constructing the following loss function to implement the blind estimation network:

[0035] Loss = ||Z1-Y1||1;

[0036] In each iteration of the network, the trainable parameters B and R are projected into the corresponding range, specifically: B∈[0,1], and R is nonnegative.

[0037] Optionally, in step 2, the estimated point spread function B and spectral response function R are obtained through training, with the learning rate set to 0.00005 and the number of iterations to 3000.

[0038] Optionally, in step 1, a low spatial resolution hyperspectral image Y and a high spatial resolution multispectral image Z are generated according to the Wald protocol.

[0039] Optionally, in step 1, the real image is spatially blurred and downsampled at a sampling rate of 4 or 8, with a corresponding spatial blur kernel size of 4*4 or 8*8 and a variance of 2, thereby generating the low spatial resolution hyperspectral image Y.

[0040] The beneficial effects of this invention are:

[0041] This application uses simulated image input to design an autoencoder network, setting B and R as learnable parameters. B and R are estimated from the simulated data, addressing the issues of scarce sensor information, difficulty in accurately adapting sensor information to data fusion models, and uncertainty and realism in degradation. Conventional methods often assume both are known, which is unrealistic. This application achieves blind hyperspectral and multispectral image fusion, making it more reasonable and realistic. Furthermore, based on matrix factorization, a deep CP decomposition module is designed to solve the problem, considering the three-dimensional information of the data, utilizing global correlation to refine image features, recover high-frequency information, and incorporating a designed regularization term to obtain richer spectral and spatial features, resulting in better fusion results. In practice, it demonstrates excellent performance. Attached Figure Description

[0042] To more clearly illustrate the technical solutions in the embodiments of the present invention, the accompanying drawings used in the description of the embodiments will be briefly introduced below. Obviously, the accompanying drawings described below are only some embodiments of the present invention. For those skilled in the art, other drawings can be obtained based on these drawings without creative effort.

[0043] Figure 1 This is a framework diagram of the unsupervised hyperspectral image super-resolution method based on matrix factorization network provided by the present invention.

[0044] Figure 2(a) is a structural diagram of the Deep CP Decomposition Module (DCPD) provided in this application, and (b) and (c) correspond to the Low Rank Tensor Generation Module (LRTG) and the Low Rank Tensor Learning Module (LRTL), respectively.

[0045] Figure 3These are residual images of the PU dataset fused according to embodiments of the present invention, wherein (a) is the label, (b) is the CNMF fused image, (c) is the HySure fused image, and (d) is the fused image of the present invention.

[0046] Figure 4 These are residual images of the CAVE dataset fused according to embodiments of this application, wherein (a) is a label, (b) is a CNMF fused image, (c) is a HySure fused image, and (d) is a fused image of this invention. Detailed Implementation

[0047] To make the objectives, technical solutions, and advantages of the present invention clearer, the embodiments of the present invention will be described in further detail below with reference to the accompanying drawings.

[0048] First, the basic theoretical knowledge involved in this invention will be introduced:

[0049] 1. Hyperspectral data often exhibits correlation between bands. A spectral vector of size Lh typically exists in a subspace smaller than Lh, and can therefore be represented as Z = EX, where E ∈ R. Lh×Ls Let Ls be a matrix whose Ls column spans the same subspace as the Z column, and X∈R. Ls×nm This is the coefficient matrix. Here, E is called the terminal membership matrix, and X is called the abundance matrix;

[0050] 2. Hyperspectral image super-resolution is performed under simulated conditions. The label data is used as a spatial response through Gaussian blur, i.e., the point spread function, to generate low-resolution data. Different sensors correspond to different spectral response functions, which are used to generate corresponding multispectral data.

[0051] Example 1:

[0052] This embodiment provides a hyperspectral image super-resolution method, including:

[0053] Step 1: Process the hyperspectral image to generate a simulated low spatial resolution hyperspectral image Y and a high spatial resolution multispectral image Z as training data;

[0054] Step 2: Input the data pair (Y, Z) generated in Step 1 into the autoencoder network, and train iteratively to obtain the point spread function B and the spectral response function R;

[0055] Step 3: For the target image, the high spatial resolution hyperspectral image X, it is assumed to be a linear combination of the terminal member matrix and the corresponding abundance matrix; the initial terminal member matrix A0 and the initial abundance matrix S0 are used as trainable parameters, and the estimated B and R in Step 2 are input into the deep CP decomposition network to solve for the desired terminal member matrix A through the deep CP decomposition network.

[0056] Step 4: Combine the point spread function B and spectral response function R obtained in Step 2, and combine the degradation model and matrix factorization model to model the target, design regularization terms, and iteratively solve the terminal member matrix A and the corresponding abundance matrix S through network expansion to finally obtain the fused image.

[0057] Example 2:

[0058] This embodiment provides an unsupervised hyperspectral image super-resolution method based on matrix factorization networks, which generally consists of four steps, as follows: Figure 1 As shown, the process includes: generating training samples; estimating the point spread function B and spectral response function R through an autoencoder network; learning the terminal member matrix A through a low-rank prior learning network; modeling the super-resolution problem by combining a degradation model and a matrix factorization model; designing regularization terms using the estimated point spread function B and spectral response function R; solving the network iteratively; and finally obtaining the terminal member matrix A and the abundance matrix S, thus obtaining the fusion result X.

[0059] Step 1: Generate a low spatial resolution hyperspectral image (LRHS) Y and a high spatial resolution multispectral image (HRMS) Z according to the Wald protocol, as training data, specifically including:

[0060] The real image is spatially blurred and downsampled in all directions with a sampling rate ratio of 4 or 8. The corresponding spatial blur kernel size is 4*4 or 8*8, and the variance is set to 2, generating a low spatial resolution hyperspectral image Y.

[0061] For different datasets, the spectral responses corresponding to different satellites were selected. The CAVE dataset used the spectral response corresponding to the D700 sensor, and the Pavia University dataset used the spectral response corresponding to ROSIS, generating a high spatial resolution multispectral image Z.

[0062] In addition, Gaussian noise was added to the low spatial resolution hyperspectral image Y (signal-to-noise ratio = 30 dB) and the high spatial resolution multispectral image (signal-to-noise ratio = 40 dB).

[0063] Step 2: Use the data pair (Y, Z) generated in Step 1 as input to the autoencoder network, and use the point spread function B and spectral response function R as trainable parameters of the network. The estimated point spread function B and spectral response function R are obtained through training. The learning rate is set to 0.00005, and the number of iterations is 3000. The specific process is as follows:

[0064] Step 2.1: Since LRHS is the spatial domain degradation of HRHS, and HRMS is the spectral degradation of HRHS, the degradation models for Y and Z can be modeled as follows:

[0065] Z = RX + N m Y = XB + N h

[0066] Where, N m and N h These represent different signal-to-noise ratios, and X represents a high spatial resolution hyperspectral image of the target.

[0067] By combining the two equations above, we can model the problem and obtain ZB = RY. Therefore, we can use an optimization algorithm to solve for min||ZB-RY|| and obtain the estimated point spread function B and spectral response function R through iteration.

[0068] Step 2.2: Solve using an autoencoder network:

[0069] Input the generated data pair (Y, Z), randomly initialize the point spread function B and the spectral response function R. The width and height of B are both (2*ratio)-1, where ratio is the sampling rate, i.e. the ratio of the sizes of the two data spaces. The size of R is (C, c), where C represents the number of bands in the hyperspectral data and c represents the number of bands in the multispectral data. Use B and R as trainable parameters of the network.

[0070] Step 2.3: Using R as the convolution kernel, tensor multiplication is implemented through convolution. This is similarly applied to ZB, ultimately optimizing the objective function min||ZB-RY||. Z1 and Y1 are iteratively solved, where Z1 and Y1 are obtained through a quadratic degradation of the learned B and R. The obtained Z1 and Y1 are normalized, and B and R are preserved. The following loss function is constructed to implement the blind estimation network:

[0071] Loss = ||Z1-Y1||1;

[0072] In addition, to achieve physical constraints, the trainable parameters B and R are projected into the corresponding ranges in each iteration of the network, specifically: B is in the interval [0,1] and R is non-negative.

[0073] Step 3: Initialize the terminal membership matrix and the corresponding abundance matrix. Use the initial terminal membership matrix A0 and the initial abundance matrix S0 as trainable parameters. Input B and R estimated in Step 2 into the network, and solve for A through the designed deep CP decomposition network. The specific steps are as follows:

[0074] Step 3.1: Based on the matrix factorization model: X = AS, establish the objective function: min||X-AS||, where A represents the terminal member matrix and S represents the corresponding abundance matrix;

[0075] Step 3.2: Randomly initialize the terminal member matrix A and the abundance matrix S to obtain the initial terminal member matrix A0 and the corresponding initial abundance matrix S0, and obtain the initial target image X0 = A0S0; where A has a size of (1, rank, H, W), where rank is the set CP rank, H and W are the length and width of the multispectral data, and S has a size of (C, rank, 1, 1), where C is the number of bands of the hyperspectral data. The point spread function B and the spectral response function R obtained in Step 2 are used as inputs.

[0076] Step 3.3: Input the initial terminal member matrix A0 into the low-rank tensor generation module (LRTG). First, extract shallow features through a convolutional Conv, and then use Conv-ReLU-Conv to generate the input feature map. The size of the spatial kernel is 3x3. Input the extracted feature map into the depth CP decomposition module. First, perform global pooling GAP on the channel, height and width dimensions respectively to obtain vectors C1, C2 and C3. Perform Kroneckor product on the three to generate tensor O. Then, introduce residual learning to learn tensors Oi (i = 1, 2... r) of different frequencies:

[0077] For the generated O i Extracting the residual can be considered as O i Unrecoverable high-frequency information. The remaining portion is processed again by the low-rank tensor generation module to generate O. i+1 (i = 1, 2, ..., r). Finally, the low-rank tensor learning module (LRTL) is entered. Given a predefined CP rank r, r tensors {O1, O2, ..., O} are generated through r iterations of the low-rank tensor generation module. r Finally, all generated tensors are aggregated into a low-rank tensor. Since higher frequency context information is more difficult to recover, this embodiment adds skip connections from the tensor that connects the beginning. This embodiment uses a 3×3 convolutional layer to aggregate the stacked tensors into the desired end member matrix A.

[0078] Step 4: Combining the point spread function B and spectral response function R obtained in Step 2, and incorporating the degradation model and matrix factorization model, the target is modeled, a regularization term is designed, and the terminal member matrix A and the corresponding abundance matrix S are iteratively solved through network expansion to finally obtain the fused image. The specific process is as follows:

[0079] From the degeneracy model Z = RX + N m Y = XB + N hModeling with the simultaneous matrix decomposition formula X = AS yields the optimization objective: min||Z-RAS||1+||Y-ASB||1, with the addition of the design regularization term ||A-A0||1. Using the point spread function B and spectral response function R estimated in step 2, along with the initialized A0 and S0, we solve for A and S:

[0080] Solve for A i With a fixed value S, this method uses a deep CP decomposition network to solve the problem. Inputs are Y, Z, and A. i-1 R, B, output A i (i = 1, 2, ..., n, where n is the iteration number, which is the number of times the network is trained).

[0081] Solve for S i By fixing A, and taking the derivative of the objective function with respect to S, we can obtain...

[0082] S i =(RA) i-1 ) -1 (Z-RA i-1 S)+A -1 (YA i-1 SB)B -1

[0083] Iterate the above steps until convergence, solve for A and S, and finally obtain the fused result X = AS. The loss function for network training is as follows, where Z2 and Y2 are the low spatial resolution hyperspectral image and high spatial resolution multispectral image obtained by further degrading the generated X. The degradation uses the real B and R, and the parameter a is set to 0.008.

[0084] Loss=||Z2-Z||1+||Y2-Y||1+a||A-A0||1

[0085] The trained network was used for validation, and the fusion quality was evaluated.

[0086] The method of the present invention will now be further described with reference to a specific embodiment. Using Pavia University (PU, see...) Figure 3 ) and CAVE (see Figure 4 Hyperspectral public datasets are fused using the method described in this invention.

[0087] The Pavia University dataset was obtained using the Reflection Optical Systems Imaging Spectrometer (ROSIS), with a spectral range of 0.43–0.86 μm. The ROSIS sensor has 115 spectral bands, and after noise removal, 103 spectral bands remain. The image size is 610 × 340 pixels, with a spatial resolution of 1.3 m per pixel. We selected the upper left 512 × 256 pixel portion as the reference image.

[0088] The CAVE dataset consists of 32 indoor HSIs captured by a high-quality generalized classification pixel camera. The HSIs have 31 bands with wavelengths ranging from 400 nm to 700 nm, and each band has a spatial size of 512 × 512.

[0089] The experimental environment used in this invention is as follows: Python 3.7, CPU: i5-8400, main frequency 2.8GHz, GPU: GTX-1060, and memory 3GB.

[0090] When the sampling rate ratio = 4, for the PU dataset, the number of endmembers is 80, the learning rate is 0.0005, and the number of iterations is 3000.

[0091] For the CAVE dataset, the number of endmembers is 12, the learning rate is 0.001, and the number of iterations is 5000.

[0092] The comparative experiment was Coupled Nonnegative Matrix Factorization Unmixing (CNMF), which evaluated the final fusion result using five metrics: Peak Signal-to-Noise Ratio (PSNR), Spectral Angular Distance (SAM), Root Mean Square Error (RMSE), Relative Global Dimensionless Error (ERGAS), and Structural Similarity (SSIM). These metrics can be used to evaluate the image reconstruction performance of the super-resolution method. Higher values ​​for PSNR and SSIM indicate better performance, while lower values ​​for SAM, ERGAS, and RMSE indicate better performance.

[0093] Table 1: Comparison of metrics for each method when the sampling rate is 4.

[0094]

[0095] The experimental results show that, for a sampling rate of 4, all indicators of the method of this invention are higher than those of CNMF and Hysure.

[0096] On the PU dataset, our method outperforms CNMF by 5.1369, 1.2655, 0.0075, 1.2827, and 0.0292 in PSNR, SAM, RMSE, ERGAS, and SSIM, respectively. Compared to Hysure, our method outperforms Hysure by 6.1967, 3.0006, 0.0075, 1.3783, and 0.0292 in PSNR, SAM, RMSE, ERGAS, and SSIM, respectively.

[0097] On the CAVE dataset, our method outperforms CNMF by 0.8394, 1.0049, 0.0007, 0.4047, and 0.0033 in PSNR, SAM, RMSE, ERGAS, and SSIM, respectively. Compared to Hysure, our method outperforms Hysure by 3.3418, 1.0049, 0.0036, 1.1702, and 0.0069 in PSNR, SAM, RMSE, ERGAS, and SSIM, respectively.

[0098] When the sampling rate ratio = 8, for the PU dataset, the number of endmembers is 80, the learning rate is 0.0005, and the number of iterations is 2400; for the CAVE dataset, the number of endmembers is 12, the learning rate is 0.0015, and the number of iterations is 10000.

[0099] Table 2: Comparison of metrics for each method when the sampling rate is 8%.

[0100]

[0101] The experimental results show that, for a sampling rate of 8, all indicators of the method of this invention are higher than those of CNMF and Hysure.

[0102] On the PU dataset, our method outperforms CNMF by 10.0993, 2.6589, 0.0260, 2.0042, and 1.47 in PSNR, SAM, RMSE, ERGAS, and SSIM, respectively. Compared to Hysure, our method outperforms Hysure by 10.1791, 5.2097, 0.0263, 1.8896, and 0.0707 in PSNR, SAM, RMSE, ERGAS, and SSIM, respectively.

[0103] On the CAVE dataset, our method outperforms CNMF by 1.0052, 0.5850, 0.0021, 0.4735, and 0.0182 in PSNR, SAM, RMSE, ERGAS, and SSIM, respectively. Compared to Hysure, our method outperforms Hysure by 1.8717, 1.6100, 0.0035, 0.5361, and 0.0388 in the same metrics.

[0104] The above results demonstrate that, under the same conditions, the method of this invention, utilizing a blind estimation network in conjunction with a deep CP decomposition module, achieves improvements in both spatial and spectral resolution of the fused image. Compared to CNMF and Hysure, this invention exhibits advantages in hyperspectral super-resolution.

[0105] Figure 4 (d) is the fused residual image when the sampling rate ratio = 4 in the Pavia University dataset. Figure 3 (d) is the residual image after fusion when the sampling rate ratio of the CAVE dataset is 4. The residual image can reflect the size of the pixel difference between the two images. The smaller the difference, the better the effect and the less obvious the residual image. It can be seen from the figure that for the CAVE dataset and the PU dataset, the effect from worst to best according to the difference of the residual image is Hysure, CNMF, and the present invention.

[0106] Some steps in the embodiments of the present invention can be implemented using software, and the corresponding software program can be stored in a readable storage medium, such as an optical disc or a hard disk.

[0107] The above description is only a preferred embodiment of the present invention and is not intended to limit the present invention. Any modifications, equivalent substitutions, improvements, etc., made within the spirit and principles of the present invention should be included within the protection scope of the present invention.

Claims

1. A hyperspectral image super-resolution method, characterized in that, The hyperspectral image super-resolution method includes: Step 1: Process the hyperspectral image to generate a simulated low spatial resolution hyperspectral image Y and a high spatial resolution multispectral image Z as training data; Step 2: Input the data pair (Y, Z) generated in Step 1 into the autoencoder network, and train iteratively to obtain the point spread function B and the spectral response function R; Step 3: For the target high spatial resolution hyperspectral image X, assume it is a linear combination of the terminal member matrix and the corresponding abundance matrix; use the initial terminal member matrix A0 and the initial abundance matrix S0 as trainable parameters, input B and R estimated in step 2 into the deep CP decomposition network, and solve for the desired terminal member matrix A through the deep CP decomposition network; Step 4: Combine the point spread function B and spectral response function R obtained in Step 2, and combine the degradation model and matrix factorization model to model the target, design regularization terms, and iteratively solve the terminal member matrix A and the corresponding abundance matrix S through network expansion to finally obtain the fused image; Step 2 specifically includes: Step 2.1: The degradation models for the low spatial resolution hyperspectral image Y and the high spatial resolution multispectral image Z are respectively modeled as follows: Z=RX+N m ,Y=XB+N h Where, N m and N h These represent different signal-to-noise ratios, and X represents a high spatial resolution hyperspectral image of the target. Combining the two equations above, we get: ZB = RY, and the final objective function is: min||ZB - RY||. Step 2.2: Solve using an autoencoder network: Input the generated data pair (Y, Z), randomly initialize the point spread function B and spectral response function R. The width and height of B are both (2*ratio)-1, where ratio is the sampling rate, i.e. the ratio of the spatial resolution of the two data points. The size of R is (C, c), where C represents the number of bands in the hyperspectral data and c represents the number of bands in the multispectral data. Use B and R as trainable parameters of the network. Step 2.3: Use R as the convolution kernel to implement tensor multiplication in the form of convolution. Similarly, apply it to ZB to optimize the objective function min||ZB-RY||. Iterate to solve Z1 and Y1, where Z1 and Y1 are obtained by the second degeneration of the learned B and R. Normalize the obtained Z1 and Y1 and save B and R. Step 3 specifically includes: Step 3.1: Based on the matrix factorization model: X=AS, establish the objective function: min||X-AS||, where A represents the terminal member matrix and S represents the corresponding abundance matrix; Step 3.2: Randomly initialize the terminal member matrix A and the abundance matrix S to obtain the initial terminal member matrix A0 and the corresponding initial abundance matrix S0, thus obtaining the initial target image X0=A0S0; where A has a size of (1, rank, H, W), where rank is the set CP rank, H and W are the length and width of the multispectral data, and S has a size of (C, rank, 1, 1), where C is the number of bands of the hyperspectral data. The point spread function B and the spectral response function R obtained in Step 2 are used as inputs. Step 3.3: Input the initial terminal member matrix A0 into the low-rank tensor generation module LRTG. First, a convolutional Conv is used to extract shallow features, and then a convolutional Conv-ReLU-Conv is used to generate the input feature map. The size of the spatial kernel is 3x3. The extracted feature map is input into the depth CP decomposition module. First, global pooling GAP is performed on the channel, height and width dimensions respectively to obtain vectors C1, C2 and C3. The three are multiplied by Kroneckor to generate tensor O. Then, residual learning is introduced to learn tensors O of different frequencies. i (i=1, 2…r): For the generated O i Extracting the residual portion, treating it as O i Unrecoverable high-frequency information; the remaining portion is processed again by the low-rank tensor generation module LRTG to generate O. i+1 (i=1, 2…r); Finally, the low-rank tensor learning module LRTL is entered. Given a predefined CP rank r, r tensors {O1, O2, ... O} are generated through r iterations of the low-rank tensor generation module. r Finally, all generated tensors are aggregated into a low-rank tensor, and convolutional layers are used to aggregate the stacked tensors into the desired end member matrix A.

2. The hyperspectral image super-resolution method according to claim 1, characterized in that, Step 4 includes: From the degeneracy model Z=RX+N m Y=XB+N h The optimization objective is obtained by modeling the simultaneous matrix decomposition formula X=AS: min ||Z-RAS||1 + ||Y-ASB||1, and the design regularization term ||A-A0||1 is added. Using the point spread function B and spectral response function R estimated in step 2, as well as the initialized A0 and S0, A and S are solved: Solve for A i Solve using a deep CP decomposition network, with input A. i-1 Output A i i = 1, 2…n, where n is the iteration number, which is the number of times the network is trained. Solve for S i S is set as a learnable parameter that is automatically learned and updated as the network is trained. Iterate through the above steps until convergence, solve for A and S, and finally obtain the fused result X=AS.

3. The hyperspectral image super-resolution method according to claim 2, characterized in that, The loss function for training the deep CP decomposition network is as follows: Loss=||Z2-Z||1 + ||Y2-Y||1+ a|| A-A0||1 Z2 and Y2 are low-spatial-resolution hyperspectral and high-spatial-resolution multispectral images obtained by further degrading the generated X. The degradation uses real B and R, and parameter a is a constant.

4. The hyperspectral image super-resolution method according to claim 2, characterized in that, The parameter a is set to 0.

008.

5. The hyperspectral image super-resolution method according to claim 1, characterized in that, Also includes: Construct the following loss function to implement the blind estimation network: Loss = ||Z1-Y1||1; In each iteration of the network, the trainable parameters B and R are projected into the corresponding range, specifically: B∈[0,1], and R is nonnegative.

6. The hyperspectral image super-resolution method according to claim 1, characterized in that, In step 2, the estimated point spread function B and spectral response function R are obtained through training, with the learning rate set to 0.00005 and the number of iterations to 3000.

7. The hyperspectral image super-resolution method according to claim 1, characterized in that, In step 1, a low spatial resolution hyperspectral image Y and a high spatial resolution multispectral image Z are generated according to the Wald protocol.

8. The hyperspectral image super-resolution method according to claim 7, characterized in that, Step 1 involves spatially blurring the real image and downsampling it at a sampling rate of 4 or 8. The corresponding spatial blur kernel size is 4*4 or 8*8, and the variance is set to 2, thereby generating the low spatial resolution hyperspectral image Y.

Images