# Single-Image Super-Resolution Using Panchromatic Gradient Prior and Variational Model

Authors: Yingying Xu; Jianhua Li; Haifeng Song; Lei Du
Journal: Mathematical Problems in Engineering
 (2021)

---

## Abstract

Single-image super-resolution (SISR) is a resolution enhancement technique and is known as an ill-posed problem. Motivated by the idea of pan-sharping, we propose a novel variational model for SISR. The structure tensor of the input low-resolution image is exploited to obtain the gradient of an imaginary panchromatic image. Then, by constraining the gradient consistency, the image edges and details can be better recovered during the procedure of restoration of high-resolution images. Besides, we resort to the nonlocal sparse and low-rank regularization of image patches to further improve the super-resolution performance. The proposed variational model is efficiently solved by ADMM-based algorithm. We do extensive experiments in natural images and remote sensing images with different magnifying factors and compare our method with three classical super-resolution methods. The subjective visual impression and quantitative evaluation indexes both show that our method can obtain higher-quality results.


---

## Body

## 1. Introduction

Image super-resolution (SR) is one of the most fundamental problems in the field of image processing, which aims to reconstruct clear and accurate high-resolution images from degraded low-resolution images. In other words, SR technology can recover a high-resolution (HR) image from one or more low-resolution (LR) input images [1]. In recent years, image SR has attracted great attention from academia and industry. The quality of the reconstructed SR image also greatly affects the accuracy of other computer vision tasks (such as image classification, image segmentation, and target detection). At present, SR technology has made great progress; it has been widely used in high-definition digital television, remote sensing monitoring, medical image reconstruction, image restoration, military reconnaissance, and other fields [2]. However, as an ill-posed problem, the existing image SR technology cannot always obtain satisfactory reconstruction results. Therefore, the study of SR is still a challenging but significant research topic.According to the number of input LR images, the SR can be classified into multi-image super-resolution (MISR) [3–5] and single-image super-resolution (SISR) [6–8]. MISR uses the complementary information provided by multiple images of the same scene to enhance resolution. However, in a practical application, it is much difficult to obtain image sequences of the same scene and achieve precise subpixel registration. Meanwhile, SISR is suitable for more scenarios, which enhances the resolution of one input image based on some prior information [9, 10]. To date, SISR methods can be broadly classified into interpolation-based methods, reconstruction-based methods, and learning-based methods [11].Interpolation methods like nearest-neighbor interpolation, bilinear interpolation, and bicubic Interpolation [12, 13] estimate the intensity at a point using the information of adjacent pixels. While these methods are easier calculated and very speedy, they tend to generate excessive smoothing and jagged artifacts. Reconstruction-based SISR methods [14–16] generally utilize image priors to restrict the possible solution space with an advantage of recovering sharp details. However, these methods are usually time-consuming, and their performance suffers a rapid degradation when the amplification factor increases.Learning-based SISR methods include example-based methods, neighbor-embedding-based methods, sparse-representation-based methods, and deep-learning-based methods [17, 18]. The basic idea of these approaches is to obtain the correspondence between the LR image and the HR image by a priori knowledge training and then reconstruct an HR image from the input LR image using the correspondence that has been learned. Freeman et al. [19] presented an example-based method with a Markov network. Chang et al. [20] proposed a neighbor-embedding-based method that needs less training set. Yang et al. [21] proposed a sparse-representation-based algorithm to obtain more favorable results and further improve their work by jointly training coupled dictionaries for LR and HR image patch pairs [22]. Very recently, deep-learning-based SISR methods, using convolutional network [23, 24], residual network [25, 26], generative adversarial network [27], attention network [28, 29], and so on, have become popular and demonstrated good performance with the rapid development of deep learning technology.Pan-sharpening [30] is a multispectral image fusion and super-resolution technique; the panchromatic (PAN) image has higher spatial resolution than the corresponding multispectral image but only has a single band. Pan-sharpening fuses PAN image and multispectral image together to enhance both spatial and spectral resolutions of the data. Motivated by the idea of pan-sharpening, we present a novel variational approach for single-image super-resolution. In this paper, we assume that there is a PAN image of the input LR image. By utilizing the structure tensor [31] of input image, we can construct the gradient of the PAN image; then we build a variational model combined with low-rank and sparse representation of similar patches to effectively fuse the constructed PAN information with LR image and obtain HR image. The proposed approach addresses the super-resolution problem of both single nature image and the multispectral data without panchromatic image.The remainder of this paper is organized as follows. In Section2, we present the proposed SISR model in detail. The numerical algorithm for solving our model is given in Section 3. The experimental results and analysis are discussed in Section 4. Finally, conclusion is presented in Section 5. In addition, the descriptions of the acronyms used in this paper are listed in Table 1.Table 1
The descriptions of the acronyms used in this paper.
AcronymDescriptionSRSuper-resolutionSISRSingle-image super-resolutionMISRMulti-image super-resolutionHRHigh resolutionLRLow resolutionPANPanchromaticADMMAlternating direction method of multipliersPSNRPeak signal-to-noise ratioRMSERoot mean square errorSSIMStructure similarity indexSCSparse-coding-based super-resolution methodSCBPIterative-back-projection-based SC method


## 2. Variational Model

Drawing on the idea of pan-sharpening, we construct a variational fusion model to realize single-image super-resolution. As panchromatic image contains the edges and details needed for resolution enhancement, we first construct the information of PAN image which does not actually exist. According to [31], the key information of the assumed PAN image can be exploited from the structure tensor of the input image. Then the gradient of PAN image can be fused into the LR image to enhance the HR details. The similar image patch pairs extracted from LR and HR images are constrained by low-rank and sparse regularization to improve the SR performance.


### 2.1. Image Degradation Model

We consider that the observed low-resolution imageH is a downsampled version of the high-resolution image F, and the degradation model [32] can be written as(1)H=DF+ɛ,where D represents a downsampling operator and ϵ is the additive Gaussian white noise.


### 2.2. Construct the Gradient of PAN Image

We construct the gradient of the PAN image corresponding to input LR image using the knowledge of structure tensor [31]. Firstly, we apply quadratic interpolation to the given LR image H to obtain the enlarged image H∗ in desired scale. Then, matrix G known as the structure tensor [31] can be expressed as(2)G=∑i=1Nsi∇Hi∗si∇Hi∗T=∑i=1Nsi∇xHi∗si∇yHi∗si∇xHi∗si∇yHi∗=∑i=1Nsi∇xHi∗2,∑i=1Nsi2∇xHi∗∇yHi∗∑i=1Nsi2∇yHi∗∇xHi∗,∑i=1Nsi∇yHi∗2,where N denotes the number of spectral bands; si denotes the weight and it can be determined by si = ∇Hi∗/∑t=1N∇Ht∗21/2.As a symmetric matrix,G can be decomposed as(3)G=QΛQT,where Q is an orthogonal matrix and the columns of which are the eigenvectors of G; Λ is a diagonal matrix and the diagonal elements of which are the eigenvalues of G. Equation (3) can also be written as(4)G=λ1θ1θ1T+λ2θ2θ2T,where λ1 is the maximum eigenvalue and it gives the maximum rate of change of H∗; λ2 is the minimum eigenvalue and it gives the minimum rate of change. The corresponding eigenvectors θ1 and θ2 give the directions of change.The PAN imageP can theoretically capture the basic geometry of the image H∗. Consider image P with λmax = ∇P2 and λmin = 0, whose structure tensor ∇P∇PT should approximately be equal to G. Thus, we have the following equation:(5)∇P∇PT=λ1θ1θ1T.In order to solve∇P, equation (5) can be rewritten as follows:(6)∇P∇PT=λ1θ1λ1θ1T.To specify the sign of the eigenvectors, the target gradient of the constructed PAN imageP is thus obtained as(7)∇P=λ1θ1signθ1⋅∑i=1Nsi∇Hi∗.


### 2.3. Variational Model Based on Gradient Consistency

In this subsection, we build a variational model to fuse the low-resolution imageH and the gradient of its corresponding PAN image P; then the model is optimized together with the nonlocal sparse and low-rank regularization to obtain the high-resolution image F.


#### 2.3.1. Consistency of Gradient

Bands of an ideal HR imageF often have similar structural information with the PAN image P; namely, bands of HR image closely approximate the PAN image in gradient information. Therefore, we can use constraint about the consistency of gradient to establish a relationship model. We use Fi to denote the i-th band of image F and then set up the gradient consistency regularization term as follows:(8)∑i=1Nαi∇Fi=∇P,where αi denotes the weight of the i-th band in F, and we generally set αi=1/N in the absence of special requirements.


#### 2.3.2. Sparse and Low-Rank Decomposition

Motivated by the effective work of optical flow estimation with nonlocal sparse and low-rank regularization in [33], we adopt a joint sparse and low-rank decomposition term formulated as follows:(9)RkF=Lk+Sk+ε,where Rk denotes the operator that extracts the k-th patch of image and groups its nonlocal similar patches together. Then the grouped similar patches are decomposed to the low-rank component Lk and the sparse component Sk, and ε denotes the Gaussian white noise.


#### 2.3.3. Total Energy

Combining the equality constraint in (1), (8), and (9), the proposed variational SISR model can be formulated as the following optimization objective:(10)minF,Lk,Sk12DF−H22+μ12∑i=1Nαi∇Fi−∇P22+∑i=1N∑k=1Mμ22RkFi−Lk−Sk22+μ32rankLk+μ42Sk0,where M is the number of image patches that are extracted from the i-th band of image F, N is the number of spectral bands, μ1,μ2,μ3,μ4 are parameters, rankLk denotes the rank of low-rank component, and ℓ0 norm denotes the sparse regularization. It is well known that the rank-solution and ℓ0 norm are discrete combinatorial optimization problems; they are both NP-hard. Therefore, in order to make the minimization tractable, we resort to their convex relaxed modification, and model (10) can be reformulated as the following convex optimization problem:(11)minF,Lk,Sk12DF−H22+μ12∑i=1Nαi∇Fi−∇P22+∑i=1N∑k=1Mμ22RkFi−Lk−Sk22+μ32Lk∗+μ42Sk1,where ⋅∗ denotes the nuclear norm, which is the convex envelope (tightest convex surrogate) of rank, and ℓ1 norm is the nearest convex norm to ℓ0 norm.


## 2.1. Image Degradation Model

We consider that the observed low-resolution imageH is a downsampled version of the high-resolution image F, and the degradation model [32] can be written as(1)H=DF+ɛ,where D represents a downsampling operator and ϵ is the additive Gaussian white noise.


## 2.2. Construct the Gradient of PAN Image

We construct the gradient of the PAN image corresponding to input LR image using the knowledge of structure tensor [31]. Firstly, we apply quadratic interpolation to the given LR image H to obtain the enlarged image H∗ in desired scale. Then, matrix G known as the structure tensor [31] can be expressed as(2)G=∑i=1Nsi∇Hi∗si∇Hi∗T=∑i=1Nsi∇xHi∗si∇yHi∗si∇xHi∗si∇yHi∗=∑i=1Nsi∇xHi∗2,∑i=1Nsi2∇xHi∗∇yHi∗∑i=1Nsi2∇yHi∗∇xHi∗,∑i=1Nsi∇yHi∗2,where N denotes the number of spectral bands; si denotes the weight and it can be determined by si = ∇Hi∗/∑t=1N∇Ht∗21/2.As a symmetric matrix,G can be decomposed as(3)G=QΛQT,where Q is an orthogonal matrix and the columns of which are the eigenvectors of G; Λ is a diagonal matrix and the diagonal elements of which are the eigenvalues of G. Equation (3) can also be written as(4)G=λ1θ1θ1T+λ2θ2θ2T,where λ1 is the maximum eigenvalue and it gives the maximum rate of change of H∗; λ2 is the minimum eigenvalue and it gives the minimum rate of change. The corresponding eigenvectors θ1 and θ2 give the directions of change.The PAN imageP can theoretically capture the basic geometry of the image H∗. Consider image P with λmax = ∇P2 and λmin = 0, whose structure tensor ∇P∇PT should approximately be equal to G. Thus, we have the following equation:(5)∇P∇PT=λ1θ1θ1T.In order to solve∇P, equation (5) can be rewritten as follows:(6)∇P∇PT=λ1θ1λ1θ1T.To specify the sign of the eigenvectors, the target gradient of the constructed PAN imageP is thus obtained as(7)∇P=λ1θ1signθ1⋅∑i=1Nsi∇Hi∗.


## 2.3. Variational Model Based on Gradient Consistency

In this subsection, we build a variational model to fuse the low-resolution imageH and the gradient of its corresponding PAN image P; then the model is optimized together with the nonlocal sparse and low-rank regularization to obtain the high-resolution image F.


### 2.3.1. Consistency of Gradient

Bands of an ideal HR imageF often have similar structural information with the PAN image P; namely, bands of HR image closely approximate the PAN image in gradient information. Therefore, we can use constraint about the consistency of gradient to establish a relationship model. We use Fi to denote the i-th band of image F and then set up the gradient consistency regularization term as follows:(8)∑i=1Nαi∇Fi=∇P,where αi denotes the weight of the i-th band in F, and we generally set αi=1/N in the absence of special requirements.


### 2.3.2. Sparse and Low-Rank Decomposition

Motivated by the effective work of optical flow estimation with nonlocal sparse and low-rank regularization in [33], we adopt a joint sparse and low-rank decomposition term formulated as follows:(9)RkF=Lk+Sk+ε,where Rk denotes the operator that extracts the k-th patch of image and groups its nonlocal similar patches together. Then the grouped similar patches are decomposed to the low-rank component Lk and the sparse component Sk, and ε denotes the Gaussian white noise.


### 2.3.3. Total Energy

Combining the equality constraint in (1), (8), and (9), the proposed variational SISR model can be formulated as the following optimization objective:(10)minF,Lk,Sk12DF−H22+μ12∑i=1Nαi∇Fi−∇P22+∑i=1N∑k=1Mμ22RkFi−Lk−Sk22+μ32rankLk+μ42Sk0,where M is the number of image patches that are extracted from the i-th band of image F, N is the number of spectral bands, μ1,μ2,μ3,μ4 are parameters, rankLk denotes the rank of low-rank component, and ℓ0 norm denotes the sparse regularization. It is well known that the rank-solution and ℓ0 norm are discrete combinatorial optimization problems; they are both NP-hard. Therefore, in order to make the minimization tractable, we resort to their convex relaxed modification, and model (10) can be reformulated as the following convex optimization problem:(11)minF,Lk,Sk12DF−H22+μ12∑i=1Nαi∇Fi−∇P22+∑i=1N∑k=1Mμ22RkFi−Lk−Sk22+μ32Lk∗+μ42Sk1,where ⋅∗ denotes the nuclear norm, which is the convex envelope (tightest convex surrogate) of rank, and ℓ1 norm is the nearest convex norm to ℓ0 norm.


## 2.3.1. Consistency of Gradient

Bands of an ideal HR imageF often have similar structural information with the PAN image P; namely, bands of HR image closely approximate the PAN image in gradient information. Therefore, we can use constraint about the consistency of gradient to establish a relationship model. We use Fi to denote the i-th band of image F and then set up the gradient consistency regularization term as follows:(8)∑i=1Nαi∇Fi=∇P,where αi denotes the weight of the i-th band in F, and we generally set αi=1/N in the absence of special requirements.


## 2.3.2. Sparse and Low-Rank Decomposition

Motivated by the effective work of optical flow estimation with nonlocal sparse and low-rank regularization in [33], we adopt a joint sparse and low-rank decomposition term formulated as follows:(9)RkF=Lk+Sk+ε,where Rk denotes the operator that extracts the k-th patch of image and groups its nonlocal similar patches together. Then the grouped similar patches are decomposed to the low-rank component Lk and the sparse component Sk, and ε denotes the Gaussian white noise.


## 2.3.3. Total Energy

Combining the equality constraint in (1), (8), and (9), the proposed variational SISR model can be formulated as the following optimization objective:(10)minF,Lk,Sk12DF−H22+μ12∑i=1Nαi∇Fi−∇P22+∑i=1N∑k=1Mμ22RkFi−Lk−Sk22+μ32rankLk+μ42Sk0,where M is the number of image patches that are extracted from the i-th band of image F, N is the number of spectral bands, μ1,μ2,μ3,μ4 are parameters, rankLk denotes the rank of low-rank component, and ℓ0 norm denotes the sparse regularization. It is well known that the rank-solution and ℓ0 norm are discrete combinatorial optimization problems; they are both NP-hard. Therefore, in order to make the minimization tractable, we resort to their convex relaxed modification, and model (10) can be reformulated as the following convex optimization problem:(11)minF,Lk,Sk12DF−H22+μ12∑i=1Nαi∇Fi−∇P22+∑i=1N∑k=1Mμ22RkFi−Lk−Sk22+μ32Lk∗+μ42Sk1,where ⋅∗ denotes the nuclear norm, which is the convex envelope (tightest convex surrogate) of rank, and ℓ1 norm is the nearest convex norm to ℓ0 norm.


## 3. Numerical Algorithm

In this section, the numerical procedure of the proposed model will be implemented. We adopt the alternating direction method of multipliers (ADMM) [34] to solve the minimization problem (11). The flexibility of the ADMM lies in the fact that it splits the initial optimization problem into several subproblems.


### 3.1. Solution of the Sparse MatrixSk

For each band, we fixF and the low-rank matrix Lk; the optimization problem of sparse matrix Sk is as follows:(12)S^k=argminSkμ22RkF−Lk−Sk22+μ42Sk1.We can obtain the solution of equation (12) by soft threshold algorithm as(13)S^k=softRkF−Lk,μ4μ2,where softa,τ=maxa22−τ,0signa.


### 3.2. Solution of the Low-Rank MatrixLk

For each band, we fixF and the sparse matrix Sk; the optimization problem of low-rank matrix Lk is as follows:(14)L^k=argminLkμ22RkF−Lk−Sk22+μ32Lk∗.We can obtain the solution of equation (14) by singular value threshold algorithm as(15)L^k=UDτ^ΣVT,Dτ^Σ=diagσi−τ^+,where σi−τ^+=maxσi−τ^,0, σi are singular values of matrix RkF−Sk, U and V are left and right singular vectors of σi, and τ^=μ3/μ2.


### 3.3. Solution of the SR ImageF

In order to separately solve the fidelity term and two regularization terms in (11), we introduce two intermediate variables U and V which are close to F. Thus, we need to iteratively solve the three following subproblems to obtain the solution of SR image.(1)
We use the intermediate variableU to replace F in 1/2DF−H22. The optimization problem for U is as follows:(16)U^=argminU12DU−H22+W1,F−U+ρ2F−U22,where⋅,⋅ denotes the inner product, W1 denotes the Lagrange multiplier matrix, and ρ is a positive penalty parameter. We set the derivative of (16) with respect to U to 0:(17)DTDU−H−W1−ρF−U=0,which yields(18)U=DTH+W1+ρFDTD+ρ.The Lagrange multiplierW1 is updated as follows:(19)W1=W1+ρF−U.(2)
We introduce an intermediate variableVi to replace Fi in the gradient consistency term. The optimization problem for Vi is as follows:(20)Vi^=argminViμ12∑i=1Nαi∇Vi−∇P22+∑i=1NW2i,Fi−Vi+ρ2∑i=1NFi−Vi22,whereW2 denotes the Lagrange multiplier matrix. We set the derivative of (20) with respect to Vi to 0:(21)μ1αi∇T∑i=1Nαi∇Vi−∇P−W2i−ρFi−Vi=0.Using the fast Fourier transform (FFT), we can obtain the closed-form solution ofVi from (21) directly:(22)Vi=ℱ−1ℱ−μ1αi∇T∑j≠iNαj∇Vj−∇P+W2i+ρFiℱμ1αi2∇T∇+ρ,whereℱ and ℱ−1 denote the FFT and inverse FFT, respectively. The Lagrange multiplier W2 is updated as follows:(23)W2=W2+ρF−V.(3)
AfterU and V are substituted into equation (11), the optimization subproblem of F becomes(24)F^=argminF∑i=1N∑k=1Mμ22RkFi−Lk−Sk22+∑i=1NW1i,Fi−Ui+∑i=1NW2i,Fi−Vi+ρ2∑i=1NFi−Ui22+ρ2∑i=1NFi−Vi22.We set the derivative of (24) with respect to Fi to 0, which yields(25)Fi=μ2RkTLk+Sk+ρUi+Vi−W1i−W2iμ2RkTRk+2ρ.Overall, taking all above analyses into account, we can summarize the complete numerical procedure for the proposed method. The detailed descriptions are shown in Algorithm1.Algorithm 1: Algorithm 1A variational approach for single-image SR.
(i)
Input:The LR imageH;Magnification factor.Output:The SR imageF∗.Initialization:Compute an initial SR imageF0 via bicubicinterpolation;Calculate the target gradient∇P via equation (7).For: t=1,2,…,TFor: i=1,2,…,N1. For each extracted image patchRkF0i, searchfor a set of similar patches across the wholeimage by patch matching.2. Decompose the grouped patches into low-rankand sparse components.(a) Calculate the sparse matrixSk via (13);(b) Calculate the low-rank matrixLk via (15).3. Reconstruct the SR imageF.(a) Update the imageF via (25);(b) Update the intermediate variableU and V via (18) and (22);(c) Update the Lagrange multiplierW1 and W2 via (19) and (23);(d)F0  =  F;End;End;F∗=F0


## 3.1. Solution of the Sparse MatrixSk

For each band, we fixF and the low-rank matrix Lk; the optimization problem of sparse matrix Sk is as follows:(12)S^k=argminSkμ22RkF−Lk−Sk22+μ42Sk1.We can obtain the solution of equation (12) by soft threshold algorithm as(13)S^k=softRkF−Lk,μ4μ2,where softa,τ=maxa22−τ,0signa.


## 3.2. Solution of the Low-Rank MatrixLk

For each band, we fixF and the sparse matrix Sk; the optimization problem of low-rank matrix Lk is as follows:(14)L^k=argminLkμ22RkF−Lk−Sk22+μ32Lk∗.We can obtain the solution of equation (14) by singular value threshold algorithm as(15)L^k=UDτ^ΣVT,Dτ^Σ=diagσi−τ^+,where σi−τ^+=maxσi−τ^,0, σi are singular values of matrix RkF−Sk, U and V are left and right singular vectors of σi, and τ^=μ3/μ2.


## 3.3. Solution of the SR ImageF

In order to separately solve the fidelity term and two regularization terms in (11), we introduce two intermediate variables U and V which are close to F. Thus, we need to iteratively solve the three following subproblems to obtain the solution of SR image.(1)
We use the intermediate variableU to replace F in 1/2DF−H22. The optimization problem for U is as follows:(16)U^=argminU12DU−H22+W1,F−U+ρ2F−U22,where⋅,⋅ denotes the inner product, W1 denotes the Lagrange multiplier matrix, and ρ is a positive penalty parameter. We set the derivative of (16) with respect to U to 0:(17)DTDU−H−W1−ρF−U=0,which yields(18)U=DTH+W1+ρFDTD+ρ.The Lagrange multiplierW1 is updated as follows:(19)W1=W1+ρF−U.(2)
We introduce an intermediate variableVi to replace Fi in the gradient consistency term. The optimization problem for Vi is as follows:(20)Vi^=argminViμ12∑i=1Nαi∇Vi−∇P22+∑i=1NW2i,Fi−Vi+ρ2∑i=1NFi−Vi22,whereW2 denotes the Lagrange multiplier matrix. We set the derivative of (20) with respect to Vi to 0:(21)μ1αi∇T∑i=1Nαi∇Vi−∇P−W2i−ρFi−Vi=0.Using the fast Fourier transform (FFT), we can obtain the closed-form solution ofVi from (21) directly:(22)Vi=ℱ−1ℱ−μ1αi∇T∑j≠iNαj∇Vj−∇P+W2i+ρFiℱμ1αi2∇T∇+ρ,whereℱ and ℱ−1 denote the FFT and inverse FFT, respectively. The Lagrange multiplier W2 is updated as follows:(23)W2=W2+ρF−V.(3)
AfterU and V are substituted into equation (11), the optimization subproblem of F becomes(24)F^=argminF∑i=1N∑k=1Mμ22RkFi−Lk−Sk22+∑i=1NW1i,Fi−Ui+∑i=1NW2i,Fi−Vi+ρ2∑i=1NFi−Ui22+ρ2∑i=1NFi−Vi22.We set the derivative of (24) with respect to Fi to 0, which yields(25)Fi=μ2RkTLk+Sk+ρUi+Vi−W1i−W2iμ2RkTRk+2ρ.Overall, taking all above analyses into account, we can summarize the complete numerical procedure for the proposed method. The detailed descriptions are shown in Algorithm1.Algorithm 1: Algorithm 1A variational approach for single-image SR.
(i)
Input:The LR imageH;Magnification factor.Output:The SR imageF∗.Initialization:Compute an initial SR imageF0 via bicubicinterpolation;Calculate the target gradient∇P via equation (7).For: t=1,2,…,TFor: i=1,2,…,N1. For each extracted image patchRkF0i, searchfor a set of similar patches across the wholeimage by patch matching.2. Decompose the grouped patches into low-rankand sparse components.(a) Calculate the sparse matrixSk via (13);(b) Calculate the low-rank matrixLk via (15).3. Reconstruct the SR imageF.(a) Update the imageF via (25);(b) Update the intermediate variableU and V via (18) and (22);(c) Update the Lagrange multiplierW1 and W2 via (19) and (23);(d)F0  =  F;End;End;F∗=F0


## 4. Experimental Results and Analysis

In this section, we first evaluate the sensitivity of the proposed algorithm with respect to the main parameters, that is,μ1, μ2, μ3, and μ4. Then, to demonstrate the SR effectiveness on both natural images and remote sensing images with different magnifying factors (2 and 4), we compare our method with the Bicubic interpolation method, sparse-coding-based method (SC) [22], and its back-projection enhanced version (SCBP) [35]. We evaluate the outcome of various methods by using quantitative indexes: peak signal-to-noise ratio (PSNR), root mean square error (RMSE), and structure similarity index (SSIM) [36]; they are calculated on RGB channels for nature image and on each band for multispectral image. The higher PSNR and SSIM values and lower RMSE value represent the better SR performance.


### 4.1. Adjustment of Parameters

The parameterμ1 is the weight to the consistency of gradient, μ2 controls the low-rank and sparse regularization of the similar patches, and μ3 and μ4 control the contribution of low-rank and sparse regularization terms. We adjust the parameter values in proper ranges and demonstrate the average quantitative metrics as functions of the parameters on the experimental dataset; then the optimal parameter combination is found with extensive numerical tests. Figure 1 shows the average PSNR, RMSE, and SSIM results when the parameter μ1 varies from 0 to 0.5 with fixed μ2 = 1.6, μ3 = 0.2, and μ4=0.2. From Figure 1, we can see that when μ1=0, all quantitative metrics get the worst values, which illustrate the effectiveness of the gradient consistency constraint. Meanwhile, Figure 1 shows that our method achieves the best performance when μ1=0.4. Figure 2 shows the average PSNR, RMSE, and SSIM results under various combinations of μ2 and μ3 with fixed μ1=0.4 and μ4=0.2. The analysis of μ2 and μ3 shows that the trend of the overall performance gets better and then worse. When the maximal PSNR is selected, we have μ2 = 1.8 and μ3 = 0.3; when the maximal SSIM is selected, we have μ2 = 1.4 and μ3 = 0.1. In order to synthesize the performance of our algorithm in these two indexes, we normalize PSNR and SSIM by introducing E = 0.8∗PSNR + 0.2∗SSIM, and E gets its maximum when μ2 = 1.6 and μ3 = 0.2. Besides, we find that the performance of the proposed method is insensitive to the parameter μ4, and we simply set μ4=0.2.Figure 1
SR results under differentμ1. (a) PSNR; (b) RMSE; (c) SSIM.
(a)(b)(c)Figure 2
SR results under differentμ2 and μ3. (a) PSNR; (b) RMSE; (c) SSIM.
(a)(b)(c)


### 4.2. SR Results and Analysis

#### 4.2.1. Natural Images Super-Resolution

To evaluate the performance of the proposed super-resolution algorithm, we use bicubic interpolation and SC and SCBP algorithms for comparison. Experimental results of the images “Flower,” “Leaf,” and “Lena” with magnification factor of 2 are shown in Table2 and Figure 3, while the results with magnification factor of 4 are shown in Table 3 and Figure 4.Table 2
SR results of different methods with magnification factor of 2.
ImagesBicubicSCSCBPOursFlowerPSNR32.735632.703234.858836.5869RMSE5.88525.90724.60893.7774SSIM0.91630.89980.92950.9351LeafPSNR23.780723.282025.853628.6582RMSE16.500917.476012.99759.4108SSIM0.86500.83860.89970.9323LenaPSNR32.053331.701934.084436.1762RMSE6.36616.62895.03873.9603SSIM0.91430.89390.92890.9319Figure 3
SR results of images “Flower,” “Leaf,” and “Lena” with magnification factor of 2. Each odd line from left to right: original images; results by bicubic interpolation; results by SC method; results by SCBP method; results by our method. Each even line from left to right: the magnification of the small red rectangle in the corresponding each example.Table 3
SR results of different methods with magnification factor of 4.
ImagesBicubicSCSCBPOursFlowerPSNR26.112328.625428.440828.8676RMSE12.61609.44659.64949.1867SSIM0.74160.79220.78560.7969LeafPSNR17.263619.322619.012919.3259RMSE34.843027.568228.568927.5578SSIM0.50810.59190.59050.6295LenaPSNR25.440827.906227.498728.0027RMSE13.630210.261910.754810.1486SSIM0.73810.78330.77610.7996Figure 4
SR results of images “Flower,” “Leaf,” and “Lena” with magnification factor of 4. Each line from left to right: original images; results by bicubic interpolation; results by SC method; results by SCBP method; results by our method.As we can see from the visual comparison, our method could obtain desired results; it effectively reduces zigzag and distortion of the recovered SR images so that it better preserves the image details compared to other competitors. From the close-ups of the regions in red boxes, we can see clearly that the SR images obtained by bicubic interpolation and SC methods are more blurred than ours, and SCBP method could restore clear edges but generate some obvious artifacts in SR images. In summary, our algorithm is superior to the contrast methods in terms of reconstructing clearer SR images with sharper edges and without artifacts. The above observations can be quantitatively confirmed by the Tables2 and 3 which record the PSNR, RMSE, and SSIM values of the SR results of the three images with magnification factors of 2 and 4, respectively. The proposed method obtains the best metrics values compared with the other three methods.


#### 4.2.2. Multispectral Images Super-Resolution

The adopted multispectral images contain four bands (red, blue, green, and infrared bands), and our algorithm is executed at each band. Same as the nature image experiments, we compared the performance of our method with bicubic interpolation, SC, and SCBP. Experimental results for the multispectral datasets “Field,” “Tree,” and “Beach” with magnification factor of 2 are shown in Figure5 and Table 4, while the results for magnification factor of 4 are shown in Table 5.Figure 5
SR results of images “Field,” “Tree,” and “Beach” with magnification factor of 2. Each odd line from left to right: original images; results by bicubic interpolation; results by SC method; results by SCBP method; results by our method. Each even line from left to right: the magnification of the small red rectangle in the corresponding each example.Table 4
SR results of different methods with magnification factor of 2.
ImagesBicubicSCSCBPOursFieldPSNR24.450424.505822.793226.8614RMSE15.276415.179218.487611.5736SSIM0.74930.74750.68720.8255TreePSNR23.526424.039022.184825.9000RMSE16.991116.016619.828912.9283SSIM0.72370.72480.66640.8139BeachPSNR22.997622.915020.901824.7080RMSE18.057718.230022.985414.8300SSIM0.74270.74010.68210.8166Table 5
SR results of different methods with magnification factor of 4.
ImagesBicubicSCSCBPOursFieldPSNR19.539219.594417.886121.0251RMSE26.889426.718932.893422.6614SSIM0.45360.45370.38250.5437TreePSNR19.039819.514817.692520.5211RMSE28.480626.956032.987724.0150SSIM0.40180.41380.35100.5144BeachPSNR18.394018.242016.446319.5784RMSE30.679031.220438.390826.7681SSIM0.40490.39460.33960.5010The observation of visual comparison is close to that of natural image experiment. Bicubic interpolation and SC algorithms tend to reconstruct overly smooth images that lack detailed information and have fuzzy image quality, while the SR images obtained by SCBP method are overenhanced and have obvious artifacts. Comparatively speaking, our algorithm can better restore the detailed information of images and make a better texture effect. By comparing the numerical results in Tables4 and 5, our algorithm also gets the best results for all indicators. It is worth noting that, because of the lack of training set, the results of SCBP algorithm are worse than bicubic interpolation. In conclusion, our proposed algorithm is superior to other algorithms compared in the subjective visual impression and objective indexes.


## 4.1. Adjustment of Parameters

The parameterμ1 is the weight to the consistency of gradient, μ2 controls the low-rank and sparse regularization of the similar patches, and μ3 and μ4 control the contribution of low-rank and sparse regularization terms. We adjust the parameter values in proper ranges and demonstrate the average quantitative metrics as functions of the parameters on the experimental dataset; then the optimal parameter combination is found with extensive numerical tests. Figure 1 shows the average PSNR, RMSE, and SSIM results when the parameter μ1 varies from 0 to 0.5 with fixed μ2 = 1.6, μ3 = 0.2, and μ4=0.2. From Figure 1, we can see that when μ1=0, all quantitative metrics get the worst values, which illustrate the effectiveness of the gradient consistency constraint. Meanwhile, Figure 1 shows that our method achieves the best performance when μ1=0.4. Figure 2 shows the average PSNR, RMSE, and SSIM results under various combinations of μ2 and μ3 with fixed μ1=0.4 and μ4=0.2. The analysis of μ2 and μ3 shows that the trend of the overall performance gets better and then worse. When the maximal PSNR is selected, we have μ2 = 1.8 and μ3 = 0.3; when the maximal SSIM is selected, we have μ2 = 1.4 and μ3 = 0.1. In order to synthesize the performance of our algorithm in these two indexes, we normalize PSNR and SSIM by introducing E = 0.8∗PSNR + 0.2∗SSIM, and E gets its maximum when μ2 = 1.6 and μ3 = 0.2. Besides, we find that the performance of the proposed method is insensitive to the parameter μ4, and we simply set μ4=0.2.Figure 1
SR results under differentμ1. (a) PSNR; (b) RMSE; (c) SSIM.
(a)(b)(c)Figure 2
SR results under differentμ2 and μ3. (a) PSNR; (b) RMSE; (c) SSIM.
(a)(b)(c)


## 4.2. SR Results and Analysis

### 4.2.1. Natural Images Super-Resolution

To evaluate the performance of the proposed super-resolution algorithm, we use bicubic interpolation and SC and SCBP algorithms for comparison. Experimental results of the images “Flower,” “Leaf,” and “Lena” with magnification factor of 2 are shown in Table2 and Figure 3, while the results with magnification factor of 4 are shown in Table 3 and Figure 4.Table 2
SR results of different methods with magnification factor of 2.
ImagesBicubicSCSCBPOursFlowerPSNR32.735632.703234.858836.5869RMSE5.88525.90724.60893.7774SSIM0.91630.89980.92950.9351LeafPSNR23.780723.282025.853628.6582RMSE16.500917.476012.99759.4108SSIM0.86500.83860.89970.9323LenaPSNR32.053331.701934.084436.1762RMSE6.36616.62895.03873.9603SSIM0.91430.89390.92890.9319Figure 3
SR results of images “Flower,” “Leaf,” and “Lena” with magnification factor of 2. Each odd line from left to right: original images; results by bicubic interpolation; results by SC method; results by SCBP method; results by our method. Each even line from left to right: the magnification of the small red rectangle in the corresponding each example.Table 3
SR results of different methods with magnification factor of 4.
ImagesBicubicSCSCBPOursFlowerPSNR26.112328.625428.440828.8676RMSE12.61609.44659.64949.1867SSIM0.74160.79220.78560.7969LeafPSNR17.263619.322619.012919.3259RMSE34.843027.568228.568927.5578SSIM0.50810.59190.59050.6295LenaPSNR25.440827.906227.498728.0027RMSE13.630210.261910.754810.1486SSIM0.73810.78330.77610.7996Figure 4
SR results of images “Flower,” “Leaf,” and “Lena” with magnification factor of 4. Each line from left to right: original images; results by bicubic interpolation; results by SC method; results by SCBP method; results by our method.As we can see from the visual comparison, our method could obtain desired results; it effectively reduces zigzag and distortion of the recovered SR images so that it better preserves the image details compared to other competitors. From the close-ups of the regions in red boxes, we can see clearly that the SR images obtained by bicubic interpolation and SC methods are more blurred than ours, and SCBP method could restore clear edges but generate some obvious artifacts in SR images. In summary, our algorithm is superior to the contrast methods in terms of reconstructing clearer SR images with sharper edges and without artifacts. The above observations can be quantitatively confirmed by the Tables2 and 3 which record the PSNR, RMSE, and SSIM values of the SR results of the three images with magnification factors of 2 and 4, respectively. The proposed method obtains the best metrics values compared with the other three methods.


### 4.2.2. Multispectral Images Super-Resolution

The adopted multispectral images contain four bands (red, blue, green, and infrared bands), and our algorithm is executed at each band. Same as the nature image experiments, we compared the performance of our method with bicubic interpolation, SC, and SCBP. Experimental results for the multispectral datasets “Field,” “Tree,” and “Beach” with magnification factor of 2 are shown in Figure5 and Table 4, while the results for magnification factor of 4 are shown in Table 5.Figure 5
SR results of images “Field,” “Tree,” and “Beach” with magnification factor of 2. Each odd line from left to right: original images; results by bicubic interpolation; results by SC method; results by SCBP method; results by our method. Each even line from left to right: the magnification of the small red rectangle in the corresponding each example.Table 4
SR results of different methods with magnification factor of 2.
ImagesBicubicSCSCBPOursFieldPSNR24.450424.505822.793226.8614RMSE15.276415.179218.487611.5736SSIM0.74930.74750.68720.8255TreePSNR23.526424.039022.184825.9000RMSE16.991116.016619.828912.9283SSIM0.72370.72480.66640.8139BeachPSNR22.997622.915020.901824.7080RMSE18.057718.230022.985414.8300SSIM0.74270.74010.68210.8166Table 5
SR results of different methods with magnification factor of 4.
ImagesBicubicSCSCBPOursFieldPSNR19.539219.594417.886121.0251RMSE26.889426.718932.893422.6614SSIM0.45360.45370.38250.5437TreePSNR19.039819.514817.692520.5211RMSE28.480626.956032.987724.0150SSIM0.40180.41380.35100.5144BeachPSNR18.394018.242016.446319.5784RMSE30.679031.220438.390826.7681SSIM0.40490.39460.33960.5010The observation of visual comparison is close to that of natural image experiment. Bicubic interpolation and SC algorithms tend to reconstruct overly smooth images that lack detailed information and have fuzzy image quality, while the SR images obtained by SCBP method are overenhanced and have obvious artifacts. Comparatively speaking, our algorithm can better restore the detailed information of images and make a better texture effect. By comparing the numerical results in Tables4 and 5, our algorithm also gets the best results for all indicators. It is worth noting that, because of the lack of training set, the results of SCBP algorithm are worse than bicubic interpolation. In conclusion, our proposed algorithm is superior to other algorithms compared in the subjective visual impression and objective indexes.


## 4.2.1. Natural Images Super-Resolution

To evaluate the performance of the proposed super-resolution algorithm, we use bicubic interpolation and SC and SCBP algorithms for comparison. Experimental results of the images “Flower,” “Leaf,” and “Lena” with magnification factor of 2 are shown in Table2 and Figure 3, while the results with magnification factor of 4 are shown in Table 3 and Figure 4.Table 2
SR results of different methods with magnification factor of 2.
ImagesBicubicSCSCBPOursFlowerPSNR32.735632.703234.858836.5869RMSE5.88525.90724.60893.7774SSIM0.91630.89980.92950.9351LeafPSNR23.780723.282025.853628.6582RMSE16.500917.476012.99759.4108SSIM0.86500.83860.89970.9323LenaPSNR32.053331.701934.084436.1762RMSE6.36616.62895.03873.9603SSIM0.91430.89390.92890.9319Figure 3
SR results of images “Flower,” “Leaf,” and “Lena” with magnification factor of 2. Each odd line from left to right: original images; results by bicubic interpolation; results by SC method; results by SCBP method; results by our method. Each even line from left to right: the magnification of the small red rectangle in the corresponding each example.Table 3
SR results of different methods with magnification factor of 4.
ImagesBicubicSCSCBPOursFlowerPSNR26.112328.625428.440828.8676RMSE12.61609.44659.64949.1867SSIM0.74160.79220.78560.7969LeafPSNR17.263619.322619.012919.3259RMSE34.843027.568228.568927.5578SSIM0.50810.59190.59050.6295LenaPSNR25.440827.906227.498728.0027RMSE13.630210.261910.754810.1486SSIM0.73810.78330.77610.7996Figure 4
SR results of images “Flower,” “Leaf,” and “Lena” with magnification factor of 4. Each line from left to right: original images; results by bicubic interpolation; results by SC method; results by SCBP method; results by our method.As we can see from the visual comparison, our method could obtain desired results; it effectively reduces zigzag and distortion of the recovered SR images so that it better preserves the image details compared to other competitors. From the close-ups of the regions in red boxes, we can see clearly that the SR images obtained by bicubic interpolation and SC methods are more blurred than ours, and SCBP method could restore clear edges but generate some obvious artifacts in SR images. In summary, our algorithm is superior to the contrast methods in terms of reconstructing clearer SR images with sharper edges and without artifacts. The above observations can be quantitatively confirmed by the Tables2 and 3 which record the PSNR, RMSE, and SSIM values of the SR results of the three images with magnification factors of 2 and 4, respectively. The proposed method obtains the best metrics values compared with the other three methods.


## 4.2.2. Multispectral Images Super-Resolution

The adopted multispectral images contain four bands (red, blue, green, and infrared bands), and our algorithm is executed at each band. Same as the nature image experiments, we compared the performance of our method with bicubic interpolation, SC, and SCBP. Experimental results for the multispectral datasets “Field,” “Tree,” and “Beach” with magnification factor of 2 are shown in Figure5 and Table 4, while the results for magnification factor of 4 are shown in Table 5.Figure 5
SR results of images “Field,” “Tree,” and “Beach” with magnification factor of 2. Each odd line from left to right: original images; results by bicubic interpolation; results by SC method; results by SCBP method; results by our method. Each even line from left to right: the magnification of the small red rectangle in the corresponding each example.Table 4
SR results of different methods with magnification factor of 2.
ImagesBicubicSCSCBPOursFieldPSNR24.450424.505822.793226.8614RMSE15.276415.179218.487611.5736SSIM0.74930.74750.68720.8255TreePSNR23.526424.039022.184825.9000RMSE16.991116.016619.828912.9283SSIM0.72370.72480.66640.8139BeachPSNR22.997622.915020.901824.7080RMSE18.057718.230022.985414.8300SSIM0.74270.74010.68210.8166Table 5
SR results of different methods with magnification factor of 4.
ImagesBicubicSCSCBPOursFieldPSNR19.539219.594417.886121.0251RMSE26.889426.718932.893422.6614SSIM0.45360.45370.38250.5437TreePSNR19.039819.514817.692520.5211RMSE28.480626.956032.987724.0150SSIM0.40180.41380.35100.5144BeachPSNR18.394018.242016.446319.5784RMSE30.679031.220438.390826.7681SSIM0.40490.39460.33960.5010The observation of visual comparison is close to that of natural image experiment. Bicubic interpolation and SC algorithms tend to reconstruct overly smooth images that lack detailed information and have fuzzy image quality, while the SR images obtained by SCBP method are overenhanced and have obvious artifacts. Comparatively speaking, our algorithm can better restore the detailed information of images and make a better texture effect. By comparing the numerical results in Tables4 and 5, our algorithm also gets the best results for all indicators. It is worth noting that, because of the lack of training set, the results of SCBP algorithm are worse than bicubic interpolation. In conclusion, our proposed algorithm is superior to other algorithms compared in the subjective visual impression and objective indexes.


## 5. Conclusion

This paper presented an efficient variational method for single-image super-resolution. Motivated by image fusion and pan-sharpening, we exploit the gradient of imaginary PAN image from the input image via structure tensor and then construct a gradient consistency constraint which is supposed to fuse HR edges and details into the optimization target image. Meanwhile, we adopt the regularization of sparse and low-rank decomposition for similar image patches to further improve the SR accuracy. The variables of our variational model are iteratively optimized by ADMM-based algorithm. Extensive experiment results demonstrate the edge recovery and artifacts suppression capabilities of our method.



---

*Source: 9944385-2021-05-03.xml*