HYPERSPECTRAL IMAGE DISTRIBUTED RESTORATION METHOD AND SYSTEM BASED ON GRAPH SIGNAL PROCESSING AND SUPERPIXEL SEGMENTATION

Abstract

Provide is a novel mixed-noise removal method for HSI with large size. First, the underlying structure of the HSI is modeled by a two-layer architecture graph. The upper layer, called a skeleton graph, is a rough graph constructed by using the modified k-nearest-neighborhood algorithm and its nodes correspond to a series of superpixels formed by HSI segmentation, which can efficiently characterize the inter-correlations between superpixels, while preserving the boundary information and reducing the computational complexity. The lower layer, called detailed graph, consists of a series of local graphs which are constructed to model the similarities between pixels. Second, based on the two-layer graph architecture, the HSI restoration problem is formulated as a series of optimization problems each of which resides on a subgraph. Third, a novel distributed algorithm is tailored for the restoration problem, by using the information interaction between the nodes of skeleton graph and subgraphs.

Claims

1. A hyperspectral image distributed restoration method based on graph signal processing and superpixel segmentation, comprising: step 1, constructing an input signal model, to linearly normalize a hyperspectral image to obtain a normalized hyperspectral image; step 2, selecting an image corresponding to a target-band from the normalized hyperspectral image, and pre-denoising the image to obtain a pre-denoised image; step 3, segmenting the pre-denoised image into a plurality of superpixel regions using simple linear iterative cluster (SLIC) superpixel segmentation algorithm; step 4, constructing a skeleton graph based on the plurality of superpixel regions; step 5, constructing a local graph for all pixels of each of the plurality of superpixel regions, to obtain a plurality of local graphs respectively corresponding to the plurality of superpixel regions; step 6, constructing sub-graphs based on the skeleton graph and the plurality of local graphs; step 7, applying an exchange distributed mode to each of the sub-graphs, to perform distributed denoising on the hyperspectral image; step 8, establishing an optimization model of each of the sub-graphs, based on the each of the sub-graphs; and step 9, repeating step 7 and step 8 through a plurality of iterations, wherein each of the plurality of iterations adopts the exchange distributed mode in the step 7 to iteratively solve the optimization model in the step 8; and in a situation that each of the sub-graphs being meeting an iterative convergence condition, the plurality of iterations are stopped, and a corresponding sub-graph of the sub-graphs meeting the iterative convergence condition is no longer updated while other sub-graphs of the sub-graphs not meeting the iterative convergence condition continue to participate the plurality of iterations.

2. The hyperspectral image distributed restoration method based on graph signal processing and superpixel segmentation according to claim 1, wherein in the step 4, the constructing the skeleton graph based on the plurality of superpixel regions, comprises: modeling the plurality of superpixel regions as a plurality of skeleton graph nodes respectively, wherein a graph signal of each skeleton graph node of the plurality of skeleton graph nodes is an average value values of all pixels in a superpixel region of the plurality of superpixel regions corresponding to the skeleton graph node, a coordinate of the skeleton graph node is a coordinate of a center pixel of the all pixels whose intensity is closest to the graph signal; and constructing the skeleton graph :=({circumflex over (V)}.sub.s,,) based on the plurality of skeleton graph nodes using a modified k-Nearest Neighbor (KNN) algorithm, where {circumflex over (V)}.sub.s represents a node set composed of the plurality of skeleton graph nodes, represents an edge set composed of edges of the skeleton graph, and represents a skeleton graph adjacency matrix, which is used to describe a correlation weight matrix of the edges of the skeleton graph, wherein the skeleton graph adjacency matrix is determined by: determining an initial skeleton graph adjacency matrix based on an element a.sub.i,j of the initial skeleton graph adjacency matrix, the element a.sub.i,j is expressed as a formula (1): $\begin{array}{cc}{a}_{i,j}=\exp \left(-\frac{{.Math.I_i-I_j.Math.}_2^2}{{}_1^2}\right)\exp \left(-\frac{{.Math.x_i-x_j.Math.}_2^2}{{}_x^2}\right),& \left(1\right)\end{array}$ where I.sub.i represents a coordinate of the i-th superpixel of all pixels of each of the plurality of superpixel regions; I.sub.j represents a coordinate of the j-th superpixel of all pixels of each of the plurality of superpixel regions; x.sub.i represents a graph signal of the i-th superpixel of all pixels of each of the plurality of superpixel regions; x.sub.j represents a graph signal of the j-th superpixel of all pixels of each of the plurality of superpixel regions; and .sub.1 and .sub.x represent two constant parameters; and setting first k relative greater values of each row of the initial skeleton graph adjacency matrix to 1, setting the rest of the row of the initial skeleton graph adjacency matrix to 0, and symmetrizing and modifying the initial skeleton graph adjacency matrix after setting into an unweighted matrix as the skeleton graph adjacency matrix .

3. The hyperspectral image distributed restoration method based on graph signal processing and superpixel segmentation according to claim 1, wherein in the step 5, the m-th superpixel region of the plurality of superpixel regions has K.sub.m pixels, and the constructing the local graph for all pixels of each of the plurality of superpixel regions, to obtain the plurality of local graphs respectively corresponding to the plurality of superpixel regions, comprises: modeling all pixels of each of the plurality of superpixel regions as pixel nodes, to thereby construct the plurality of local graphs :=(,,), wherein represents a node set composed of all pixel nodes in the m-th superpixel region, represents an edge set composed of edges for connecting each pixel node of the all pixel nodes in the m-th superpixel region to surrounding four neighbors of the pixel node, and represents an unweighted four-neighborhood adjacency matrix; wherein a graph signal of each of the plurality of local graphs is a vector x.sub.m=[x.sub.m,1, x.sub.m,2, . . . , x.sub.m,K.sub.m].sup.T .sup.K.sup.m, a graph signal of an image corresponding to the b-th band of the hyperspectral image is X.sub.m,b=[X.sub.m(1,b),X.sub.m(2,b), . . . ,X.sub.m(K.sub.m,b)].sup.% .sup.K.sup.m, b=1, . . . , p, and the local graph is applied to an image corresponding to each band of the hyperspectral image; wherein a matrix is formed based on graph signals of all bands, as a graph signal of the hyperspectral image, expressed as: X.sub.m=[X.sub.m(i,b)]=[X.sub.m,1, . . . , X.sub.m,b, . . . , X.sub.m,p)].sup.K.sup.m.sup.p, where X.sub.m(i,b) represents the i-th graph signal of the image corresponding to the b-th band of the hyperspectral image; X.sub.m,b represents the m-th graph signal of the image corresponding to the b-th band of the hyperspectral image, and .sup.K.sup.m.sup.p represents a dimension K.sub.mp of the matrix.

4. The hyperspectral image distributed restoration method based on graph signal processing and superpixel segmentation according to claim 3, wherein in the step 6, a hyperspectral restoration problem is transformed into solving of a plurality of small-scale problems, each sub-graph of the plurality of small-scale problems is located on one of the sub-graphs, and each sub-graph is centered on one of the plurality of skeleton graph nodes of the skeleton graph; wherein the f-th sub-graph of the sub-graphs is expressed as a formula (2):
:=(V.sub.f,E.sub.f), f=1, 2, . . . , M.sub.s(2) where $\stackrel{\_}{V_f}={\left\{\right\}}_{m\left(f,d\right)},$ represents the f-th sub-graph, V.sub.f represents a node set composed of all sub-graph nodes of the f-th sub-graph, E.sub.f represents an edge set composed of edges between the all sub-graph nodes of the f-th sub-graph, M.sub.s represents of a total number of the sub-graphs, and f represents an integer; wherein (f, d):={j{circumflex over (V)}.sub.s:(f, j)d} represents a set of neighbor sub-graph nodes of a sub-graph node f of the f-th sub-graph within d hops; wherein (f, j) represents a shortest path from the sub-graph node f to a sub-graph node j; wherein each sub-graph comprises at least two of the plurality of local graphs, the m-th local graph of the at least two of the plurality of local graphs is one sub-graph node of the set of neighbor sub-graph nodes of the sub-graph node f corresponding to the f-th sub-graph within the d hops, that is m(f,d); wherein d represents the set of neighbor sub-graph nodes within the d hops, m represents the m-th local graph and is a sub-graph node of the f-th sub-graph, represents the node set composed of all pixel nodes in the m-th superpixel region of the plurality of superpixel regions, (f,d) represents the set of neighbor sub-graph nodes of the sub-graph node f corresponding to the f-th sub-graph of the skeleton graph within the d hops; wherein the edge set E.sub.f is obtained by: connecting a boundary node to four neighborhood nodes of the boundary node, and forming an edge between the pixel nodes of a local graph and the local graph , if there is a boundary connection for the local graph and the local graph in corresponding superpixel regions thereof; or, forming only the edge between the pixel nodes of the local graph and the local graph , if there is no boundary connection for the local graph and the local graph in corresponding superpixel regions thereof, that is, there is edge connection between the pixel nodes of the local graph and the local graph ; wherein a Laplacian matrix L.sub.f of the f-th sub-graph is obtained, and a same sub-graph is applied to all bands of the hyperspectral image corresponding to the sub-graph; wherein a graph signal of the sub-graph node V.sub.f of the f-th sub-graph is mapped into a vector X.sub.f,b=[X.sub.f(1,b),X.sub.f(2,b), . . . , X.sub.f(K.sub.m, b(].sup.T .sup.K.sup.m, where m(f, d), b=1, . . . , p, X.sub.f(i,b) represents a graph signal of the i-th sub-graph node corresponding to the b-th band in the f-th sub-graph , X.sub.f(1,b) represents a graph signal of a first sub-graph node corresponding to the b-th band in the f-th sub-graph , K.sub.m represents a total number of pixels in the superpixel region corresponding to the m-th local graph on the f-th sub-graph , and b represents the b-th band; and wherein a matrix is formed based on graph signals of all bands using a formula expressed as follows: X.sub.f=[X.sub.f(i,b)]=[X.sub.f,1, . . . , X.sub.f,b, . . . , X.sub.f,p)].sup.(K.sup.m.sup.)p, m(f,d), and thus the graph signal of the hyperspectral image on the sub-graph , f=1, 2, . . . , M.sub.s is obtained.

5. The hyperspectral image distributed restoration method based on graph signal processing and superpixel segmentation according to claim 1, wherein in the step 7, a plurality of computing units are set, each of the plurality of computing units corresponds to one of the plurality of local graphs, and a transmission link with neighbored computing units of the computing unit is established by the computing unit based on edges of the skeleton graph, and wherein the applying the exchange distributed mode to each of the sub-graphs to perform distributed denoising on the hyperspectral image, comprises: receiving computing results from the neighbored computing units when each of plurality of iterations is performed; fusing and averaging, by the computing unit, the received computing results from the neighbored computing units with a last calculation result of the computing unit; determining, by the computing unit, an optimized computing result of the computing unit and optimized computing results of the neighbored computing units based on the optimization model of each of the sub-graphs; and transmitting, by the computing unit, the optimized computing result of the computing unit and a corresponding one of the optimized computing results of the neighbored computing units to a corresponding one of the neighbored computing units.

6. The hyperspectral image distributed restoration method based on graph signal processing and superpixel segmentation according to claim 1, wherein in the step 8, the optimization model is expressed as a formula (3): $\begin{array}{cc}\begin{array}{c}\min \\ X_{f,S_f}\left(K_m\right)p\end{array}{.Math.X_f.Math.}_{*}+{.Math.X_f.Math.}_{\mathrm{GLR}}+{.Math.S_f.Math.}_{1}s.t.{.Math.O_f-X_f-S_f.Math.}_F^2,\mathrm{rank}\left(S_f\right)r,& \left(3\right)\end{array}$ where m(f,d), K.sub.m represents a total number of pixels in the m-th superpixel region, and are two parameters, and used to adjust proportions of a kernel norm, GLR and an L.sub.1 norm, S.sub.f represents a model of sparse noise of in the f-th sub-graph, O.sub.f represents hyperspectral image data originally input in the f-th sub-graph, represents a threshold value, r is a constant for indicating constraint on a rank of a hyperspectral signal; f represents an integer; wherein GLR of signals of one band of the HSI is expressed as a formula (4):
x.sub.GLR=x.sup.Tx=.sub.ij.sub.i,j.Math.(x.sub.ix.sub.j).sup.2(4), the GLR of signals of all bands of the HSI is expressed as a formula (5):
X.sub.GLR=Tr(X.sup.TX)=x.sub.1.sub.GLR+ . . . +x.sub.b.sub.GLR+ . . . +x.sub.p.sub.GLR(5), where represents a Laplacian matrix of the f-th sub-graph, .sub.i,j represents the (i,j)-th element of an adjacency matrix of the f-th sub-graph; x.sub.i represents the i-th graph signal of the f-th sub-graph; x.sub.b.sub.GLR represents a GLR corresponding to the b-th band; and Tr(.Math.) represents a trace of a matrix.

7. The hyperspectral image distributed restoration method based on graph signal processing and superpixel segmentation according to claim 6, wherein in the step 8, the optimization model is obtained according to augmented Lagrange method; wherein the formula (3) is equivalently converted into a formula (6): $\begin{array}{cc}\begin{array}{c}\min \\ H_f,X_{f,S_f}\left(K_m\right)p\end{array}{.Math.H_f.Math.}_{*}+{.Math.X_f.Math.}_{\mathrm{GLR}}+{.Math.S_f.Math.}_{2}s.t.{.Math.O_f-H_f-S_f.Math.}_F^2,\mathrm{rank}\left(H_f\right)r,H_f=X_f,& \left(6\right)\end{array}$ where H.sub.f represents an introduced auxiliary variable, and H.sub.f=X.sub.f; wherein an augmented Lagrange function of the formula (6) is expressed as a formula (7): $\left(H_f,X_f,S_f,C_f,D_f\right)={.Math.H_f.Math.}_{*}+{.Math.X_f.Math.}_{\mathrm{GLR}}+{.Math.S_f.Math.}_2+.Math.C_f,O_f-H_f-S_f.Math.+.Math.D_f,X_f-H_f.Math.+\frac{}{2}\left({.Math.O_f-H_f-S_f.Math.}_F^2+{.Math.X_f-H_f.Math.}_F^2\right)$ $s.t.\mathrm{rank}\left(H_f\right)r,$ where r represents a penalty variable, which is greater than zero, and C.sub.f and D.sub.f represent Lagrange multipliers; wherein the formula (7) is alternately optimized through following updatings: 1) updating H.sub.f $\begin{array}{cc}{H}_f^{\left(k+1\right)}=\begin{array}{c}\mathrm{argmin}\\ \mathrm{rank}\left(H_f\right)r\end{array}\left(H_f,X_f^{\left(k\right)},S_f^{\left(k\right)},C_f^{\left(k\right)},D_f^{\left(k\right)}\right)=\begin{array}{c}\mathrm{argmin}\\ \mathrm{rank}\left(H_f\right)r\end{array}{.Math.H_f.Math.}_{*}+.Math.C_f^{\left(k\right)},O_f-H_f-S_f^{\left(k\right)}.Math.+.Math.D_f^{\left(k\right)},X_f^{\left(k\right)}-H_f.Math.+\frac{}{2}\left({.Math.O_f-H_f-S_f^{\left(k\right)}.Math.}_F^2\frac{}{2}{.Math.X_f^{\left(k\right)}-H_f.Math.}_F^2\right)=\begin{array}{c}\mathrm{argmin}\\ \mathrm{rank}\left(H_f\right)r\end{array}{.Math.H_f.Math.}_{*}+({.Math.H_f-\frac{1}{2}\left(O_f+X_f^{\left(k\right)}-S_f^{\left(k\right)}+\frac{1}{}\left(C_f^{\left(k\right)}+D_f^{\left(k\right)}\right)\right).Math.}_F^2,& \left(8\right)\end{array}$ where (k+1) represents a (k+1)-th iteration, X.sub.f.sup.(k) represents the denoised hyperspectral image after the k-th iteration of the f-th sub-graph, S.sub.f.sup.(k) represents sparse noise after the k-th iteration of the f-th sub-graph, C.sub.f.sup.(k) and D.sub.f.sup.(k) represent Lagrange multipliers after the k-th iteration of the f-th sub-graph; wherein for a given matrix Z .sup.Kp with a rank r, singular value decomposition is used to decompose the given matrix, to obtain a formula (9):
Z=UE.sub.rV.sup.H,E.sub.r=diag(.sub.i)(9), where U .sup.Kr and V .sup.pr are orthogonal matrices, E.sub.r represents a diagonal matrix composed of singular values, a singular value .sub.i(1ir)0, a threshold value being greater than or equal to zero is given, and a singular value contraction operator is defined as a formula (10):
D.sub.(Z)=UD.sub.(E.sub.r)V.sup.H(10)
D.sub.(E.sub.r)=diag(max(,0)) where =.sub.i, represents a vector composed of a difference between the i-th singular value .sub.i and a threshold value , .sub.i represents the i-th singular value after singular value decomposition, represents a constant, which is the threshold value, D.sub.(.Math.) represents a solution of a singular value contraction operator, the singular value contraction operator (10) is a solution of a problem (11): $\begin{array}{cc}{D}_{}\left(Z\right)=\begin{array}{c}\mathrm{argmin}\\ \mathrm{rank}\left(H\right)r\end{array}{.Math.H.Math.}_{*}+\frac{1}{2}{.Math.H-Z.Math.}_F^2,& \left(11\right)\end{array}$ wherein the singular value contraction operator (10) is used to optimize a solution of the formula (8) to obtain a formula (12): $\begin{array}{cc}{H}_f^{\left(k+1\right)}=D_{\frac{1}{2}}\left(\frac{1}{2}\left(O_f+X_f^{\left(k\right)}-S_f^{\left(k\right)}+\frac{1}{}\left(C_f^{\left(k\right)}+D_f^{\left(k\right)}\right)\right)\right),& \left(12\right)\end{array}$ 2) updating X.sub.f $\begin{array}{cc}\begin{array}{c}{X}_f^{\left(k+1\right)}=\begin{array}{c}\mathrm{argmin}\\ X\end{array}\left(H_f^{\left(k+1\right)},X_f,S_f^{\left(k\right)},C_f^{\left(k\right)},D_f^{\left(k\right)}\right)\\ =\begin{array}{c}\mathrm{argmin}\\ X_f\end{array}{.Math.X_f.Math.}_{\mathrm{GLR}}+.Math.D_f^{\left(k\right)},X_f-H_f^{\left(k+1\right)}.Math.+\frac{}{2}{.Math.X_f-H_f^{\left(k+1\right)}.Math.}_F^2\\ =\begin{array}{c}\mathrm{argmin}\\ X_f\end{array}{.Math.X_f.Math.}_{\mathrm{GLR}}+\frac{}{2}{.Math.X_f-H_f^{\left(k+1\right)}+\frac{D_f^{\left(k\right)}}{}.Math.}_F^2,\end{array},& \left(13\right)\end{array}$ where B.sub.f=H.sub.f.sup.(k+1)D.sub.b.sup.(k)/+[b.sub.1, . . . , b.sub.p], according to the matrix form of GLR in the formula (5), the formula (13) is equivalent to a formula (14): $\begin{array}{cc}{X}_f^{\left(k+1\right)}=\begin{array}{c}\mathrm{argmin}\\ X_f=\left[x_1,.Math.,x_p\right]\end{array}{.Math.}_{j=1}^p({.Math.x_j.Math.}_{\mathrm{GLR}}+\frac{}{2}\left({.Math.x_j-b_j.Math.}_2^2\right),& \left(14\right)\end{array}$ wherein the formula (14) is decomposed into p sub-problems and solved band by band, and the j-th sub-problem of the p sub-problems is expressed as a formula (15): $\begin{array}{cc}{x}_j^{\left(k+1\right)}=\begin{array}{c}\mathrm{argmin}\\ x_j\end{array}{.Math.x_j.Math.}_{\mathrm{GLR}}+\frac{}{2}({.Math.x_j-b_j.Math.}_2^2=x_j^T{L}_{\stackrel{\_}{{}_f}}{x}_j+\frac{}{2}{.Math.x_j-b_j.Math.}_2^2,& \left(15\right)\end{array}$ wherein the formula (15) is a least squares problem, and a solution of the formula (15) is expressed as a formula (16):
x.sub.j.sup.(k+1)=(2+I).sup.1b.sub.j(16), wherein a solution of the formula (14) is obtained by combining the solutions (16) of each band and is expressed as a formula (17): $\begin{array}{cc}\left[x_1^{\left(k+1\right)},.Math.,x_f^{\left(k+1\right)},.Math.,x_p^{\left(k+1\right)}\right]=\left[{\left(2L_{\stackrel{\_}{{}_f}}+I\right)}^{-1}{b}_1,.Math.{m\left(2L_{\stackrel{\_}{{}_f}}+I\right)}^{-1}{b}_f,.Math.,{\left(2L_{\stackrel{\_}{{}_f}}+I\right)}^{-1}{b}_p\right]={\left(2L_{\stackrel{\_}{{}_f}}+I\right)}^{-1}\left[b_1,.Math.,b_f,.Math.,b_p\right],& \left(17\right)\end{array}$ wherein the formula (17) is rewritten in matrix form expressed in a formula (18):
X.sub.f.sup.(k+1)=(2+I).sup.1B.sub.f(18), wherein by derivation, the formula (13) is equivalently transformed to a formula (19): $\begin{array}{cc}\begin{array}{c}{X}_f^{\left(k+1\right)}=\begin{array}{c}\mathrm{argmin}\\ X_f\end{array}{.Math.X_f.Math.}_{\mathrm{GLR}}+\frac{}{2}{.Math.X_f-H_f^{\left(k+1\right)}+\frac{D_f^{\left(k\right)}}{}.Math.}_F^2\\ =\begin{array}{c}\mathrm{argmin}\\ X_f\end{array}\mathrm{Tr}\left(X_f^T{L}_{\stackrel{\_}{{}_f}}{X}_f\right)+\frac{}{2}{.Math.X_f-H_f^{\left(k+1\right)}+\frac{D_f^{\left(k\right)}}{}.Math.}_F^2\\ ={\left(2L_{\stackrel{\_}{{}_f}}+I\right)}^{-1}\left(H_f^{\left(k+1\right)}-\frac{D_f^{\left(k\right)}}{}\right),\end{array}& \left(19\right)\end{array}$ 3) updating S.sub.f $\begin{array}{cc}\begin{array}{c}{S}_f^{\left(k+1\right)}=\begin{array}{c}\mathrm{argmin}\\ S_f\end{array}\left(H_f^{\left(k+1\right)},X_f^{\left(k+1\right)},S_f,C_f^{\left(k\right)},D_f^{\left(k\right)}\right)\\ =\begin{array}{c}\mathrm{argmin}\\ S_f\end{array}{.Math.S_f.Math.}_1+.Math.C_f^{\left(k\right)},O_f-H_f^{\left(k+1\right)}-S_f.Math.+\frac{}{2}{.Math.O_f-H_f^{\left(k+1\right)}-S_f.Math.}_F^2\\ =\begin{array}{c}\mathrm{argmin}\\ S_f\end{array}{.Math.S_f.Math.}_1+\frac{}{2}{.Math.S_f-\left(O_f-H_f^{\left(k+1\right)}+\frac{C_f^{\left(k\right)}}{}\right).Math.}_F^2,\end{array}& \left(20\right)\end{array}$ wherein a soft threshold operator is introduced and expressed as a formula (21): $\begin{array}{cc}{R}_{}\left(\right)=\{\begin{array}{c}-,\mathrm{if}>\\ +,\mathrm{if}<\\ 0,\mathrm{otherwise}\end{array},& \left(21\right)\end{array}$ where , >0, and a solution of the formula (20) is expressed as a formula (22): $\begin{array}{cc}{S}_f^{\left(k+1\right)}=R_{\frac{}{}}(\left(O_f-H_f^{\left(k+1\right)}+\frac{C_f^{\left(k\right)}}{}\right),& \left(22\right)\end{array}$ 4) updating the Lagrange multipliers C.sub.f and D.sub.f, and the penalty variable
C.sub.f.sup.(k+1)=C.sub.f.sup.(k)+(O.sub.fH.sub.f.sup.(k+1)S.sub.f.sup.(k+1))
D.sub.f.sup.(k+1)=D.sub.f.sup.(k)+(K.sub.f.sup.(k+1)H.sub.f.sup.(k+1))(23)
.sup.(k+1)=min(.sup.(k)+.sub.max), where represents a growth step, and .sub.max represents a constant.

8. The hyperspectral image distributed restoration method based on graph signal processing and superpixel segmentation according to claim 1, wherein in the step 9, the iterative convergence condition is determined according to a formula (24):
O.sub.fH.sub.f.sup.(k+1)S.sub.f.sup.(k+1).sub.F/O.sub.f.sub.F.sub.1(24)
H.sub.f.sup.(k+1)X.sub.f.sup.(k+1).sub..sub.2, where .sub.1 and .sub.2 represent two constants, H.sub.f.sup.(k+1) represents introduced auxiliary variable and is a solution of H.sub.f=X.sub.f after (k+1) iterations, S.sub.f.sup.(k+1) represents a solution of sparse noise after (k+1) iterations of the f-th sub-graph, and X.sub.f.sup.(k+1) represents a solution of a clean hyperspectral image after (k+1) iterations of the f-th sub-graph.

9. A hyperspectral image distributed restoration system based on graph signal processing and superpixel segmentation, comprising: a memory, a processor, and a computer program stored in the memory and executable by the processor, wherein the hyperspectral image distributed restoration method based on graph signal processing and superpixel segmentation according to claim 1 is implemented when the computer program is executed by the processor.

Description

BRIEF DESCRIPTION OF THE DRAWINGS

[0085] In order to make the objectives, the technical solutions and the beneficial effects of the present disclosure clearer, the present disclosure provides the following drawings for illustration.

[0086] FIG. 1 illustrates a clean original image.

[0087] FIG. 2 illustrates an image generated after a noise is added to the clean original image.

[0088] FIG. 3 illustrates an image generated after low-rank matrix recovery (LRMR) denoising is performed on the clean original image.

[0089] FIG. 4 illustrates an image generated after total-variation-regularized low-rank matrix factorization (LRTV) denoising is performed on the clean original image.

[0090] FIG. 5 illustrates an image generated after SSSB denoising is performed on the clean original image.

[0091] FIG. 6 illustrates an image generated after hyperspectral image distributed restoration method based on graph signal processing and superpixel segmentation (GLRSSDA) denoising is performed according to an embodiment of the present disclosure.

[0092] FIG. 7 illustrates a comparison view of peak signal to noise ratios (PSNR) of each band according to an embodiment of the present disclosure.

[0093] FIG. 8 illustrates a comparison view of structural similarities (SSIM) of each band according to an embodiment of the present disclosure.

DETAILED DESCRIPTION OF EMBODIMENTS

[0094] The present disclosure is further described hereinafter combined with accompanying drawings and specific embodiments, so that those skilled in the art can better understand the present disclosure and implement the present disclosure, but the specific embodiments hereinafter are not intended to limit the present disclosure.

First Embodiment

[0095] The first embodiment of the present disclosure provides a large-scale HSI distributed restoration method based on graph signal processing and superpixel segmentation. The method may include steps 1 through step 9 as follows.

[0096] In step 1, an input signal model is constructed, to linearly normalize an HSI to a range from 0 to 1, to obtain a normalized HSI. The constructed input signal model Y=X+S+N, where Y presents an input noise signal; X represents a clean original image; S represents sparse noise, which may be an impulse noise, a bad dot, a stripe noise; and N represents a Gaussian noise.

[0097] In step 2, an image corresponding to a target-band (for example, a specified band or a random band) from the normalized HSI is selected, and the image is pre-denoised using a median filter to obtain a pre-denoised image.

[0098] In step 3, the pre-denoised image is segmented into multiple superpixel regions using simple linear iterative cluster (SLIC) superpixel segmentation algorithm.

[0099] In step 4, a skeleton graph is constructed based on the multiple superpixel regions. It is assumed that there are M.sub.s superpixel regions after the segmenting. M.sub.s superpixel regions are modeled as M.sub.s skeleton graph nodes respectively. A graph signal of each skeleton graph node is an average value of all pixels in a superpixel region of the multiple superpixel regions corresponding to the skeleton graph node, a coordinate of the skeleton graph node is a coordinate of a target pixel in the superpixel region, and an intensity value of the target pixel is closest to the graph signal. The skeleton graph :=({circumflex over (V)}.sub.s,,) is constructed based on the multiple skeleton graph nodes using a modified k-Nearest Neighbor (KNN) algorithm, where {circumflex over (V)}.sub.s represents a node set composed of the multiple skeleton graph nodes, represents an edge set composed of edges of the skeleton graph, and represents a skeleton graph adjacency matrix, which is used to describe a correlation weight matrix of the edges of the skeleton graph, and the skeleton graph adjacency matrix is determined by:

[0100] determining an initial skeleton graph adjacency matrix based on an element a.sub.i,j of the initial skeleton graph adjacency matrix, the element a.sub.i,j is expressed as a formula (1):

[00017]$\begin{array}{cc}{a}_{i,j}=\exp \left(-\frac{{.Math.I_i-I_j.Math.}_2^2}{{}_1^2}\right)\exp \left(-\frac{{.Math.x_i-x_j.Math.}_2^2}{{}_x^2}\right),& \left(1\right)\end{array}$

[0101] where I.sub.i represents a coordinate of the i-th superpixel of all pixels of each of the multiple superpixel regions; I.sub.j represents a coordinate of the j-th superpixel of all pixels of each of the multiple superpixel regions; x.sub.i represents a graph signal of the i-th superpixel of all pixels of each of the multiple superpixel regions; x.sub.j represents a graph signal of the j-th superpixel of all pixels of each of the multiple superpixel regions; and .sub.1 and .sub.x represent two constant parameters; and

[0102] setting first k relative greater values of each row of the initial skeleton graph adjacency matrix to 1, setting the rest of the row of the initial skeleton graph adjacency matrix to 0, and symmetrizing and modifying the initial skeleton graph adjacency matrix after setting into an unweighted matrix as the skeleton graph adjacency matrix .

[0103] In step 5, a local graph is constructed for all pixels of each of the multiple superpixel regions, to obtain multiple local graphs respectively corresponding to the multiple superpixel regions. The m-th superpixel region of the multiple superpixel regions has K.sub.m pixels. All pixels of each of the multiple superpixel regions are modeled as pixel nodes, to thereby construct the multiple local graphs :=(,,). represents a node set composed of all pixel nodes in the m-th superpixel region, represents an edge set composed of edges for connecting each pixel node of the all pixel nodes in the m-th superpixel region to surrounding four neighbors of the pixel node, and represents an unweighted four-neighborhood adjacency matrix, which is a sparse matrix. A graph signal of each of the multiple local graphs is a vector x.sub.m=[x.sub.m,1,x.sub.m,2, . . . ,x.sub.m,K.sub.m].sup.T .sup.K.sup.m, a graph signal of an image corresponding to the b-th band of the hyperspectral image is X.sub.m,b=[X.sub.m(1,b), X.sub.m(2,b), . . . , X.sub.m(K.sub.m, b)].sup.T .sup.K.sup.m, b=1, . . . , p, and the local graph is applied to an image corresponding to each band of the hyperspectral image. A matrix is formed based on graph signals of all bands, as a graph signal of the hyperspectral image, expressed as: X.sub.m=[X.sub.m(i,b)]=[X.sub.m,1, . . . , X.sub.m,b, . . . , X.sub.m,p)].sup.K.sup.m.sup.p.

[0104] X.sub.m(i,b) represents the i-th graph signal of the image corresponding to the b-th band of the hyperspectral image; X.sub.m,b represents the m-th graph signal of the image corresponding to the b-th band of the hyperspectral image, and .sup.K.sup.m.sup.p represents a dimension K.sub.mp of the matrix.

[0105] In step 6, sub-graphs are constructed based on the skeleton graph and the multiple local graphs. Since each sub-graph corresponds to one computing node, a hyperspectral restoration problem is transformed into solving of multiple small-scale problems, each sub-graph of the multiple small-scale problems is located on one of the sub-graphs, and each sub-graph is centered on one of the multiple skeleton graph nodes of the skeleton graph. The f-th sub-graph of the sub-graphs is expressed as a formula (2):

:=(V.sub.f,E.sub.f,), f=1, 2, . . . , M.sub.s(2).

[00018]$\stackrel{\_}{V_f}={\left\{\right\}}_{m\left(f,d\right)},$

represents the f-th sub-graph, V.sub.f represents a node set

[0106] composed of all sub-graph nodes of the f-th sub-graph, E.sub.f represents an edge set composed of edges between the all sub-graph nodes of the f-th sub-graph, M.sub.s represents of a total number of the sub-graphs, and f represents an integer.

[0107] (f, d):={j{circumflex over (V)}.sub.s:(f, j)d} represents a set of neighbor sub-graph nodes of a sub-graph node f of the f-th sub-graph within d hops; (f, j) represents a shortest path (i.e., a quantity of the hopes) from the sub-graph node f to a sub-graph node j. Each sub-graph includes multiple local graphs, the m-th local graph of the at least two of the multiple local graphs is one sub-graph node of the set of neighbor sub-graph nodes of the sub-graph node f corresponding to the f-th sub-graph within the d hops, that is m(f, d); where d represents the set of neighbor sub-graph nodes within the d hops, m represents the m-th local graph and is a sub-graph node of the f-th sub-graph, represents the node set composed of all pixel nodes in the m-th superpixel region of the multiple superpixel regions, (f,d) represents the set of neighbor sub-graph nodes of the sub-graph node f corresponding to the f-th sub-graph of the skeleton graph within the d hops.

[0108] The edge set E.sub.f is obtained by: connecting a boundary node to four neighborhood nodes of the boundary node, and forming an edge between the pixel nodes of a local graph and the local graph , if there is a boundary connection for the local graph and the local graph in corresponding superpixel regions thereof; or, forming only the edge between the pixel nodes of the local graph and the local graph , if there is no boundary connection for the local graph and the local graph in corresponding superpixel regions thereof, that is, there is edge connection between the pixel nodes of the local graph and the local graph .

[0109] A Laplacian matrix L.sub.f of the f-th sub-graph is obtained, and a same sub-graph is applied to all bands of the hyperspectral image corresponding to the sub-graph.

[0110] A graph signal of the sub-graph node V.sub.f of the f-th sub-graph is mapped into a vector X.sub.f,b=[X.sub.f(1,b),X.sub.f(2,b), . . . , X.sub.f(K.sub.m,b)].sup.T .sup.K.sup.m, where m(f,d),b=1, . . . , p, X.sub.f(i,b) represents a graph signal of the i-th sub-graph node corresponding to the b-th band in the f-th sub-graph , X.sub.f(1,b) represents a graph signal of a first sub-graph node corresponding to the b-th band in the f-th sub-graph , K.sub.m represents a total number of pixels in the superpixel region corresponding to the m-th local graph on the f-th sub-graph , and b represents the b-th band.

[0111] A matrix is formed based on graph signals of all bands using a formula expressed as follows: X.sub.f=[X.sub.f(i,b)]=[X.sub.f,1, . . . , X.sub.f,b, . . . , X.sub.f,p)].sup.(K.sup.m.sup.)p, m(f,d), and thus the graph signal of the hyperspectral image on the sub-graph , f=1, 2, . . . , M.sub.s is obtained.

[0112] In step 7, a receiving-transmitting-fusing (also referred to as exchange) distributed mode is applied to each of the sub-graphs, to perform distributed denoising on the HSI. It is assumed that there are M.sub.s computing units, each computing unit corresponds to one of the multiple local graphs , and a transmission link with neighbored computing units of the computing unit is established by the computing unit based on edges of the skeleton graph .

[0113] M.sub.s computing units perform following steps during each iteration of the optimization model of each of the sub-graphs: receiving computing results from the neighbored computing units when each of multiple iterations is performed; fusing and averaging, by the computing unit, the received computing results from the neighbored computing units with a last calculation result of the computing unit; determining, by the computing unit, an optimized computing result of the computing unit and optimized computing results of the neighbored computing units based on the optimization model of each of the sub-graphs; and transmitting, by the computing unit, the optimized computing result of the computing unit and a corresponding one of the optimized computing results of the neighbored computing units to a corresponding one of the neighbored computing units. The above steps are repeated for iteration solving of the optimization model of each of the sub-graphs.

[0114] In step 8, an optimization model of each of the sub-graphs is established based on the each of the sub-graphs; the optimization model is expressed as a formula (3):

[00019]$\begin{array}{cc}\underset{X_{f,S_f}\left(K_m\right)p}{\min }{.Math.X_f.Math.}_{*}+{.Math.X_f.Math.}_{\mathrm{GLR}}+{.Math.S_f.Math.}_{1}s.t.{.Math.O_f-X_f-S_f.Math.}_F^2,\mathrm{rank}\left(X_f\right)r.& \left(3\right)\end{array}$

[0115] m(f, d), K.sub.m represents a total number of pixels in the m-th superpixel region, and are two parameters, and used to adjust proportions of a kernel norm, GLR and an L.sub.1 norm, S.sub.f represents a model of sparse noise of in the f-th sub-graph, O.sub.f represents hyperspectral image data originally input in the f-th sub-graph, represents a threshold value, r is a constant for indicating constraint on a rank of a hyperspectral signal; f represents an integer.

[0116] GLR of signals of one band of the HSI is expressed as a formula (4):

x.sub.GLR=x.sup.Tx=.sub.ij.sub.i,j.Math.(x.sub.ix.sub.j).sup.2(4).

[0117] The GLR of signals of all bands of the HSI is expressed as a formula (5):

X.sub.GLR=Tr(X.sup.TX)=x.sub.1.sub.GLR+ . . . +x.sub.b.sub.GLR+ . . . +x.sub.p.sub.GLR(5).

[0118] Therefore, it can be regarded as the combination vector form of GLR. GLR is a quadratic form of matrix, which has many excellent characteristics. A solution of quadratic problem can be solved by simple standard linear algebra. In addition, because of the rank constraint in the optimization model, the formula (3) is non-convex, and the augmented Lagrange multiplier method (ALM) can be used to solve this problem.

[0119] represents a Laplacian matrix of the f-th sub-graph, .sub.i,j represents the (i, j)-th element of an adjacency matrix of the f-th sub-graph; x.sub.i represents the i-th graph signal of the f-th sub-graph; x.sub.b.sub.GLR represents a GLR corresponding to the b-th band; and Tr(.Math.) represents a trace of a matrix.

[0120] The formula (3) is equivalently converted into a formula (6):

[00020]$\begin{array}{cc}\underset{H_{f,},X_{f,S_f}\left(K_m\right)p}{\min }{.Math.H_f.Math.}_{*}+{.Math.X_f.Math.}_{\mathrm{GLR}}+{.Math.S_f.Math.}_{1}s.t.{.Math.O_f-H_f-S_f.Math.}_F^2,\mathrm{rank}\left(H_f\right)r,H_f=X_f,& \left(6\right)\end{array}$

[0121] where H.sub.f represents an introduced auxiliary variable, and H.sub.f=X.sub.f.

[0122] An augmented Lagrange function of the formula (6) is expressed as a formula (7):

[00021]$\begin{array}{cc}\left(H_f,X_f,S_f,C_f,D_f\right)={.Math.H_f.Math.}_{*}+{.Math.X_f.Math.}_{\mathrm{GLR}}+{.Math.S_f.Math.}_1+.Math.C_f,O_f-H_f-S_f.Math.+.Math.D_f,X_f-H_f.Math.+\frac{}{2}\left({.Math.O_f-H_f-S_f.Math.}_F^2+{.Math.X_f-H_f.Math.}_F^2\right)& \left(7\right)\end{array}$$s.t.\mathrm{rank}\left(H_f\right)r,$

[0123] where r represents a penalty variable, which is greater than zero, and C.sub.f and D.sub.f represent Lagrange multipliers. When the formula (7) is alternately optimized, one variable can be optimized while other variables are remain unchanged.

[0124] The formula (7) is alternately optimized through following updatings:

[0125] 1) updating H.sub.f

[00022]$\begin{array}{cc}{H}_f^{\left(k+1\right)}=\begin{array}{c}\mathrm{argmin}\\ \mathrm{rank}\left(H_f\right)r\end{array}\left(H_f,X_f^{\left(k\right)},S_f^{\left(k\right)},C_f^{\left(k\right)},D_f^{\left(k\right)}\right)=\begin{array}{c}\mathrm{argmin}\\ \mathrm{rank}\left(H_f\right)r\end{array}{.Math.H_f.Math.}_{*}+.Math.C_f^{\left(k\right)},O_f-H_f-S_f^{\left(k\right)}.Math.+.Math.D_f^{\left(k\right)},X_f^{\left(k\right)}-H_f.Math.+\frac{}{2}\left({.Math.O_f-H_f-S_f^{\left(k\right)}.Math.}_F^2+\frac{}{2}{.Math.X_f^{\left(k\right)}-H_f.Math.}_F^2\right)& \left(8\right)\end{array}$$=\begin{array}{c}\mathrm{argmin}\\ \mathrm{rank}\left(H_f\right)r\end{array}{.Math.H_f.Math.}_{*}+({.Math.H_f-\frac{1}{2}\left(O_f+X_f^{\left(k\right)}-S_f^{\left(k\right)}+\frac{1}{}\left(C_f^{\left(k\right)}+D_f^{\left(k\right)}\right)\right).Math.}_F^2,$

[0126] (k+1) represents a (k+1)-th iteration, X.sub.f.sup.(k) represents the denoised hyperspectral image after the k-th iteration of the f-th sub-graph, S.sub.f.sup.(k) represents sparse noise after the k-th iteration of the f-th sub-graph, C.sub.f.sup.(k) and D.sub.f.sup.(k) represent Lagrange multipliers after the k-th iteration of the f-th sub-graph.

[0127] For a given matrix Z .sup.Kp with a rank r, singular value decomposition is used to decompose the given matrix, to obtain a formula (9):

Z=UE.sub.rV.sup.H,E.sub.r=diag(.sub.i)(9),

[0128] U .sup.Kr and V .sup.pr are orthogonal matrices, E.sub.r represents a diagonal matrix composed of singular values, a singular value .sub.i(1ir)0, a threshold value being greater than or equal to zero is given, and a singular value threshold value (also referred to as ingular value contraction) operator is defined as a formula (10):

D.sub.(Z)=UD.sub.(E.sub.r)V.sup.H(10)

D.sub.(E.sub.r)=diag(max(,0))

[0129] =.sub.i, represents a vector composed of a difference between an i-th singular value .sub.i and a threshold value , .sub.i represents the i-th singular value after singular value decomposition, represents a constant, which is the threshold value, D.sub.(.Math.) represents a solution of a singular value contraction operator, the singular value contraction operator (10) is a solution of a problem (11):

[00023]$\begin{array}{cc}{D}_{}\left(X\right)=\begin{array}{c}\mathrm{argmin}\\ \mathrm{rank}\left(H\right)r\end{array}{.Math.H.Math.}_{*}+\frac{1}{2}{.Math.H-Z.Math.}_F^2.& \left(11\right)\end{array}$

[0130] The singular value contraction operator (10) is used to optimize a solution of the formula (8) to obtain a formula (12):

[00024]$\begin{array}{cc}{H}_f^{\left(k+1\right)}=D_{\frac{1}{2}}\left(\frac{1}{2}\left(O_f+X_f^{\left(k\right)}-S_f^{\left(k\right)}+\frac{1}{}\left(C_f^{\left(k\right)}+D_f^{\left(k\right)}\right)\right)\right),& \left(12\right)\end{array}$

[0131] 2) updating X.sub.f

[00025]$\begin{array}{cc}\begin{array}{c}{X}_f^{\left(k+1\right)}=\begin{array}{c}\mathrm{argmin}\\ X\end{array}\left(H_f^{\left(k+1\right)},X_f,S_f^{\left(k\right)},C_f^{\left(k\right)},D_f^{\left(k\right)}\right)\\ =\begin{array}{c}\mathrm{argmin}\\ X_f\end{array}{.Math.X_f.Math.}_{\mathrm{GLR}}+.Math.D_f^{\left(k\right)},X_f-H_f^{\left(k+1\right)}.Math.+\frac{}{2}{.Math.X_f-H_f^{\left(k+1\right)}.Math.}_F^2\\ =\begin{array}{c}\mathrm{argmin}\\ X_f\end{array}{.Math.X_f.Math.}_{\mathrm{GLR}}+\frac{}{2}{.Math.X_f-H_f^{\left(k+1\right)}+\frac{D_f^{\left(k\right)}}{}.Math.}_F^2.\end{array},& \left(13\right)\end{array}$

[0132] For convenience, B.sub.f=H.sub.f.sup.(k+1)D.sub.B.sup.(k)/+[b.sub.1, . . . , b.sub.p], and according to the matrix form of GLR in the formula (5), the formula (13) is equivalent to a formula (14):

[00026]$\begin{array}{cc}{X}_f^{\left(k+1\right)}=\begin{array}{c}\mathrm{argmin}\\ X_f=\left[x_1,.Math.,x_p\right]\end{array}{.Math.}_{j=1}^p({.Math.x_j.Math.}_{\mathrm{GLR}}+\frac{}{2}\left({.Math.x_j-b_j.Math.}_2^2\right),& \left(14\right)\end{array}$

[0133] Then, the formula (14) is decomposed into p sub-problems and solved band by band, and the j-th sub-problem of the p sub-problems is expressed as a formula (15):

[00027]$\begin{array}{cc}{x}_j^{\left(k+1\right)}=\begin{array}{c}\mathrm{argmin}\\ x_j\end{array}{.Math.x_j.Math.}_{\mathrm{GLR}}+\frac{}{2}({.Math.x_j-b_j.Math.}_2^2=x_j^T{L}_{\stackrel{\_}{{}_f}}{x}_j+\frac{}{2}{.Math.x_j-b_j.Math.}_2^2.& \left(15\right)\end{array}$

[0134] The formula (15) is a least squares problem, and a solution of the formula (15) is expressed as a formula (16):

x.sub.j.sup.(k+1)=(2+I).sup.1b.sub.j(16).

[0135] Therefore, a solution of the formula (14) can be obtained by combining the solutions (16) of each band and is expressed as a formula (17):

[00028]$\begin{array}{cc}\left[x_1^{\left(k+1\right)},.Math.,x_f^{\left(k+1\right)},.Math.,x_p^{\left(k+1\right)}\right]=\left[{\left(2L_{\stackrel{\_}{{}_f}}+I\right)}^{-1}{b}_1,.Math.,{\left(2L_{\stackrel{\_}{{}_f}}+I\right)}^{-1}{b}_f,.Math.,{\left(2L_{\stackrel{\_}{{}_f}}+I\right)}^{-1}{b}_p\right]={\left(2L_{\stackrel{\_}{{}_F}}+I\right)}^{-1}\left[b_1,.Math.,b_f,.Math.,b_p\right],& \left(17\right)\end{array}$

[0136] It is apparent that the formula (17) can be rewritten in matrix form expressed in a formula (18):

X.sub.f.sup.(k+1)=(2+I).sup.1B.sub.f(18)

[0137] By derivation, the formula (13) is equivalently transformed to a formula (19):

[00029]$\begin{array}{cc}\begin{array}{c}{X}_f^{\left(k+1\right)}=\begin{array}{c}\mathrm{argmin}\\ X_f\end{array}{.Math.X_f.Math.}_{\mathrm{GLR}}+\frac{}{2}{.Math.X_f-H_f^{\left(k+1\right)}+\frac{D_f^{\left(k\right)}}{}.Math.}_F^2\\ =\begin{array}{c}\mathrm{argmin}\\ X_f\end{array}\mathrm{Tr}\left(X_f^T{L}_{\stackrel{\_}{{}_f}}{X}_f\right)+\frac{}{2}{.Math.X_f-H_f^{\left(k+1\right)}+\frac{D_f^{\left(k\right)}}{}.Math.}_F^2\\ ={\left(2L_{\stackrel{\_}{{}_f}}+I\right)}^{-1}\left(H_f^{\left(k+1\right)}-\frac{D_f^{\left(k\right)}}{}\right).\end{array}& \left(19\right)\end{array}$

[0138] 3) updating S.sub.f

[00030]$\begin{array}{cc}\begin{array}{c}{S}_f^{\left(k+1\right)}=\begin{array}{c}\mathrm{argmin}\\ S_f\end{array}\left(H_f^{\left(k+1\right)},X_f^{\left(k+1\right)},S_f,C_f^{\left(k\right)},D_f^{\left(k\right)}\right)\\ =\begin{array}{c}\mathrm{argmin}\\ S_f\end{array}{.Math.S_f.Math.}_1+.Math.C_f^{\left(k\right)},O_f-H_f^{\left(k+1\right)}-S_f.Math.+\frac{}{2}{.Math.O_f-H_f^{\left(k+1\right)}-S_f.Math.}_F^2\\ =\begin{array}{c}\mathrm{argmin}\\ S_f\end{array}{.Math.S_f.Math.}_1+\frac{}{2}{.Math.S_f-\left(O_f-H_f^{\left(k+1\right)}+\frac{C_f^{\left(k\right)}}{}\right).Math.}_F^2.\end{array},& \left(20\right)\end{array}$

[0139] A soft threshold operator is introduced and expressed as a formula (21):

[00031]$\begin{array}{cc}{R}_{}\left(\right)=\{\begin{array}{c}-,\mathrm{if}>\\ +,\mathrm{if}<\\ 0,\mathrm{otherwise}\end{array},& \left(21\right)\end{array}$

[0140] where , >0, and a solution of the formula (20) is expressed as a formula (22):

[00032]$\begin{array}{cc}{S}_f^{\left(k+1\right)}=R_{\frac{}{}}(\left(O_f-H_f^{\left(k+1\right)}+\frac{C_f^{\left(k\right)}}{}\right),& \left(22\right)\end{array}$

[0141] 4) updating the Lagrange multipliers C.sub.f and D.sub.f, and the penalty variable

C.sub.f.sup.(k+1)=C.sub.f.sup.(k)+(O.sub.fH.sub.f.sup.(k+1)S.sub.f.sup.(k+1))

D.sub.f.sup.(k+1)=D.sub.f.sup.(k)+(X.sub.f.sup.(k+1)H.sub.f.sup.(k+1))(23)

.sup.(k+1)=min(.sup.(k)+.sub.max),

[0142] where represents a growth step, and .sub.max represents a constant.

[0143] In step 9, the step 7 and the step 8 are repeated through multiple iterations, where each of the multiple iterations adopts the receiving-transmitting-fusing distributed mode in the step 7 to iteratively solve the optimization model in the step 8; and in a situation that each of the sub-graphs being meeting an iterative convergence condition, the repeating is stopped, and a corresponding sub-graph of the sub-graphs meeting the iterative convergence condition is no longer updated while other sub-graphs of the sub-graphs not meeting the iterative convergence condition continue to participate the repeating. The iterative convergence condition is determined according to a formula (24):

O.sub.fH.sub.f.sup.(k+1)S.sub.f.sup.(k+1).sub.F/O.sub.f.sub.F.sub.1(2)

H.sub.f.sup.(k+1)X.sub.f.sup.(k+1).sub..sub.2,

[0144] where .sub.1 and .sub.2 represent two constants, H.sub.f.sup.(k+1) represents introduced auxiliary variable and is a solution of H.sub.f=X.sub.f after (k+1) iterations, S.sub.f.sup.(k+1) represents a solution of sparse noise after (k+1) iterations of the f-th sub-graph, and X.sub.f.sup.(k+1) represents a solution of a clean hyperspectral image after (k+1) iterations of the f-th sub-graph.

[0145] Further, a hyperspectral image distributed restoration system based on graph signal processing and superpixel segmentation is provided according to an embodiment of the present disclosure, including: a memory, a processor, and a computer program stored in the memory and executable by the processor, where the hyperspectral image distributed restoration method based on graph signal processing and superpixel segmentation described above is implemented when the computer program is executed by the processor.

Second Embodiment

[0146] In the second embodiment, the input hyperspectral image is a clean hyperspectral image of a Washington Shopping Center, namely X, which has 191 bands and a dimension thereof is 1280307191. Noises are added artificially, and the added noises include parse noise (including a salt and pepper noise and a deadline noise) and a Gaussian noise. An output after noise addition is Y. A three-dimensional (3D) hyperspectral image after noise addition is input into the model of this method of the present disclosure to obtain the denoised hyperspectral image.

[0147] In this embodiment, parameters of the specific steps are set as follows. In the step (2), a specific hyperspectral band is selected, and hyperspectral data of the 50-th band is taken; and a window adopted by the median filtering is [3, 3]. For parameters of the SLIC superpixel segmentation in the step (3), a size of the regions is set to 30, and a regularization parameter is set to 0.05. In the step (4), .sub.1=100, .sub.x=0.2, and k=3. In the step (6), d=1. In the step (8), =0.5, =0.23, =.sub.1 =.sub.2=1e6, r=3, =1.5, .sub.max=1e6, each of H.sub.f.sup.(0), X.sub.f.sup.(0), S.sub.f.sup.(0), C.sub.f.sup.(0), D.sub.f.sup.(0) is an all-zero matrix, and .sup.(0)=1e2. Finally, a peak signal-to-noise ratio (PSNR) and a structural similarity (SSIM) parameters of the denoised hyperspectral image X and the clean image without adding noise are calculated to evaluate the superiority of the denoised effect.

[0148] The denoised 3D hyperspectral image is input into the technical solution model with set parameters of the present disclosure. The solution GLRSSDA is compared with a hyperspectral image restoration method based on low-rank matrix recovery (LRMR) proposed by Hongyan Zhang, the hyperspectral image restoration method based on total-variation-regularized low-rank matrix factorization (LRTV) proposed by Wei He and the hyperspectral image restoration method based on superpixel segmentation of smooth band (SSSB) proposed by Yaru Fan, and the comparison graphs before and after denoising are shown in FIGS. 1 to 6, which are grayscale graphs showing a depth of black, the darker the color, the smaller a gray value of a gray scale image. FIG. 1 illustrates a clean original image, FIG. 2 illustrates an image generated after a noise is added to the clean original image, FIG. 3 illustrates an image denoised by the LRMR method, FIG. 4 illustrates an image denoised by the LRTV method, FIG. 5 illustrates an image denoised by the SSSB method, and FIG. 6 illustrates an image denoised by the GLRSSDA of the present disclosure. According to the comparison between FIGS. 1 and 2, it can be observed that more serious noise is added in this experiment, which can better reflect the effectiveness of the technical solution of the present disclosure. According to FIGS. 3 and 5, it can be observed that there is detail loss in the image, and there is over-smoothing in FIG. 4, while FIG. 6 shows a good denoising ability, which not only removes noise well, but also excels in edge details, that is to say, it can be observed that the model of the technical solution of the present disclosure has a good denoising effect. Comparing the PSNR and SSIM parameters of each band after denoising in the technical solution of the present disclosure with the solutions of the above-mentioned prior art, FIGS. 7 and 8 are obtained. It can be seen from the two drawings that the circled solid line in the drawings represents a curve of the GLRSSDA of the technical solution of the present disclosure, the solid line with . represents a curve of the LRMR, the dotted line with + represents a curve of the LRMR, and the dotted line with a five-pointed star represents a curve of the SSSB. There is a footnote corresponding to the line type and technical solution in an upper left of each of the two drawings. The higher the PSNR and SSIM, the better the denoising effect. From FIGS. 7 and 8, it can be observed that the denoising effect of PSNR and SSIM in the technical solution of the present disclosure is better than that of LRMR, LRTV and SS SB in most bands.

[0149] In summary, it can be observed that the technical solution of the present disclosure is superior to the technical solutions of the above-mentioned prior arts. The PSNR and SSIM of each band in the technical solution of the present disclosure and the above-mentioned existing solution are averaged to obtain MPSNR and MSSIM indexes, to obtain results shown in a table 1 below, and it can be observed that the technical solution of the present disclosure is superior to the above-mentioned existing solutions.

TABLE-US-00001 TABLE 1 Comparison table of denoising effect of different methods Comparison of Denoising method denoising effect LRMR LRTV SSSB GLRSSDA Parameter MPSNR 28.93 37.30 33.09 40.38 MSSIM 0.5131 0.7905 0.6615 0.8893

[0150] The large-scale hyperspectral image distributed restoration method based on graph signal processing and superpixel segmentation provided by the embodiment of the present disclosure overcomes the problem that hyperspectral images (HSI) are often disturbed by various noises.

[0151] The above-mentioned embodiments are merely preferred embodiments for fully explaining the present disclosure, and the scope of protection of the present disclosure is not limited thereto. Equivalent substitutions or changes made by those skilled in the art on the basis of the present disclosure are within the scope of protection of the present disclosure. The scope of protection of the disclosure is subject to the claims.