CN102842115B - Based on the compressed sensing image super-resolution rebuilding method of double dictionary study - Google Patents

Abstract

The present invention relates to the compressed sensing image super-resolution rebuilding method based on double dictionary study, it comprises the steps: redundant dictionary, encoder dictionary parameter training, autoregressive model weighting parameter is trained, with the redundant dictionary trained, encoder dictionary, autoregressive model weighting parameter, super-resolution rebuilding is carried out to single frames low-resolution image, this algorithm has rebuilds effective feature, is applicable to medical imaging, satellite remote sensing remote measurement, military surveillance and the numerous areas such as location and city security protection.

Description

Based on the compressed sensing image super-resolution rebuilding method of double dictionary study

Technical field:

The invention belongs to digital image processing techniques field, specifically a kind of image super-resolution rebuilding method, the high-fidelity for single image is amplified.

Background technology:

Along with popularizing of security device, the application of watch-dog also can be more and more extensive, in our life, role is also more and more important, due to monitoring the larger or object scene in the visual field from watch-dog distance more away from, the pixel be assigned on each scenery is fewer, causes the loss in detail of object scene, image thickens, be unfavorable for subsequent treatment and target identification, therefore, super-resolution rebuilding carried out to object scene and just seems particularly important.Super-resolution rebuilding technology is widely used in medical imaging, satellite remote sensing remote measurement, military surveillance and the field such as location and city security protection.

Mainly contain three major types reconstruction algorithm at present:

One class is traditional multi-frame processing algorithm, this type of algorithm is rebuild by the different information of the multiple image merging Same Scene, this type of algorithm is all based on same model, by choosing different bound term, an ill-conditioning problem is solved, this type of algorithm input requirements is higher, and the accuracy requirement of parameter estimation is very high, therefore applies limited.

Interpolation algorithm is that a class is the simplest, the highest ageing algorithm, comprise bilinear interpolation, bicubic interpolation, cubic spline interpolation etc., interpolation algorithm does not consider image immanent structure (edge), easily cause image blurring, this type of algorithm is generally as the pre-service of other super-resolution rebuilding algorithms.

Algorithm based on machine learning: such algorithm, first by obtaining the information needed for super-resolution rebuilding to image library learning training, then utilizes the information learnt to rebuild.Tradition can realize large multiple based on the method for learn-by-example and rebuild, but requires that the image in the image of input and image library is a class image, is generally applied to the image reconstruction that face etc. is special.Adopt the method for two redundant dictionary rarefaction representation can obtain good result when quality of input image is higher, but to fuzzyyer, the image reconstruction capabilities that noise is larger is just not enough.

The present invention adopts single redundancy dictionary to carry out rarefaction representation, and adopt the method for iterative shrinkage to carry out target super resolution image to Reconstructed equation and solve, reconstructed results is good.

Summary of the invention:

The object of this invention is to provide a kind of image super-resolution rebuilding method, it can keep the sharpness of image while enlarged image.

The technical solution used in the present invention is as follows:

One, redundant dictionary, encoder dictionary parameter training:

Definition for redundant dictionary, Ψ=[ψ ₁, ψ ₂..., ψ _n] ∈ R ^{m × n}for encoder dictionary, m, n are positive integer.Double dictionary refers to redundant dictionary and encoder dictionary.

The first step: super-resolution image in reading images storehouse, transfers super-resolution coloured image to gray level image, is then divided into size to be image fritter sample, the image fritter obtained is comply with left-to-right, from top to bottom, by row reading manner formed column vector, use s _i∈ R ⁿ, i=1,2 ..., Q represents that the column vector that each fritter is formed, Q are the number of total column vector;

Second step: calculate s _ivariance Var (s _i), only retain Var (s _i) be greater than the vector of threshold value TH, finally obtain training sample set S=[s ₁, s ₂... s _m];

3rd step: formula (1) is solved, adopt alternative manner to solve redundant dictionary Φ and encoder dictionary Ψ, if Θ is sparse coefficient, λ, η are constant, represent and ask l ₂norm, || ₁represent and ask l ₁norm:

$\{\mathrm{\Φ},\mathrm{\Ψ},\mathrm{\Θ}\}=\arg \underset{\mathrm{\Φ},\mathrm{\Ψ},\mathrm{\Θ}}{\min }\{{\left|\right|S-\mathrm{\Φ\Θ}\left|\right|}_2^2+\mathrm{\η}{\left|\right|\mathrm{\Θ}-\mathrm{\ΨS}\left|\right|}_2^2+\mathrm{\λ}{\left|\mathrm{\Θ}\right|}_1\}---\left(1\right)$

(1) with gaussian random matrix initialisation redundant dictionary Φ, with unit matrix initialization codes dictionary Ψ, with full null matrix initialization sparse coefficient Θ, if iterations k=0, maximum iteration time is Max_Iter, and iteration convergence controlling elements are ε;

(2) (T is defined _ζ[O]) _{i, j}=sign (O _{i, j}) max{|O _{i, j}|-ζ, 0} are threshold operation operator, and ζ represents threshold operation variable, and O represents threshold operation matrix variables, O _{i, j}be designated as the element of (i, j) under in representing matrix O, sign () is symbol manipulation operator, gets σ _Θ=2|| Φ ^tΦ+η I|| _f, wherein || || _frepresent and ask Frobenius norm, I is unit matrix, uses formula (2) to upgrade current Θ value;

${\mathrm{\Θ}}^{k+1}=T_{\mathrm{\λ}/2{\mathrm{\σ}}_{\mathrm{\Θ}}}[(1-\frac{\mathrm{\η}}{{\mathrm{\σ}}_{\mathrm{\Θ}}}){\mathrm{\Θ}}^k+\frac{1}{{\mathrm{\σ}}_{\mathrm{\Θ}}}({\mathrm{\Φ}}^T(S-\mathrm{\Φ}{\mathrm{\Θ}}^k)+\mathrm{\η\ΨS})]---\left(2\right)$

Wherein Θ ^k+1, Θ ^krepresent the Θ value of iteration kth+1 and kth step respectively, Φ ^trepresent the transposition of Φ.

(3) defining operation π (d)=d/max (1, || d||), d is vector, and this operation represents by vector projection to unit length, definition σ _Φ=2|| Θ Θ ^t|| _f, use formula (3) to upgrade current Φ value:

${\mathrm{\Φ}}^{k+1}=\mathrm{\π}({\mathrm{\Φ}}^k+\frac{1}{{\mathrm{\σ}}_{\mathrm{\Φ}}}(S-{\mathrm{\Φ}}^k\mathrm{\Θ}){\mathrm{\Θ}}^T)---\left(3\right)$

π () expression herein carries out unit length projection, wherein Φ to each row of Φ ^k+1, Φ ^krepresent the Φ value of iteration kth+1 and kth step respectively, Θ ^trepresent the transposition of Θ.

(4) σ is calculated _Ψ=2||SS ^t|| _f, use formula (4) to upgrade current Ψ value:

${\mathrm{\Ψ}}^{k+1}=\mathrm{\π}({\mathrm{\Ψ}}^k+\frac{1}{{\mathrm{\σ}}_{\mathrm{\Ψ}}}(\mathrm{\Θ}-{\mathrm{\Ψ}}^{k}S)S^T)---\left(4\right)$

π () expression herein carries out unit length projection to each row of Ψ, wherein Ψ ^k+1, Ψ ^krepresent the Ψ value of iteration kth+1 and kth step respectively, S ^trepresent the transposition of S.

(5) iterations k=k+1;

(6) (2) to (5) are repeated until the value arriving maximum iterations or formula (1) changes enough hour stopping iteration.Export Φ value and Ψ value.

Two, autoregressive model weighting parameter training:

The first step: super-resolution image in reading images storehouse, transfer image to gray level image, then with one low frequency Gaussian convolution core carries out convolution, obtain low-frequency image, then original image and low-frequency image are subtracted each other, the high-frequency information (here by its called after high frequency imaging) of difference reflection image, high frequency imaging being divided into size is image fritter, find and s in super-resolution image _ithe image fritter of corresponding position, and comply with left-to-right, from top to bottom, form column vector by row reading manner, be designated as all S=[s ₁, s ₂... s _m] corresponding high frequency fritter vector set is combined into

Second step: use K-means sorting algorithm by S ^hbe divided into K class { C ₁, C ₂... C _k, m _krepresent C _kthe number of middle vector, uses formula (5) to calculate the barycenter μ of each class _k, k=1,2 ..., K, according to S ^hclassification results, by S=[s ₁, s ₂... s _m] be also divided into K class, be expressed as { S ₁, S ₂... S _k}:

$\begin{array}{cc}{\mathrm{\μ}}_k=\frac{1}{m_k}\underset{i=1}{\overset{m_k}{\mathrm{\Σ}}}{s}_i^h& s_i^h\∈C_k\end{array}---\left(5\right)$

3rd step: if s ' _irepresent s _icenter pixel value, q _ifor s _imiddle s ' _ithe vector of neighborhood territory pixel value composition, use least square method calculating formula in a _k, identical process is carried out to all classification, obtains autoregressive model weighting parameter combination { a ₁, a ₂... a _k;

4th step: output from regression model weighting parameter combination { a ₁, a ₂... a _kand the barycenter { μ of each class ₁, μ ₂... μ _k;

Three, image super-resolution rebuilding

Redundant dictionary Φ, encoder dictionary Ψ, autoregressive model weighting parameter combination { a ₁, a ₂... a _k, the barycenter { μ of each class ₁, μ ₂... μ _kfor training in advance obtains, once training can use always.

The first step: read in the low-resolution image Y needing to carry out rebuilding, if gray level image, uses bi-cubic interpolation Y to be interpolated into the size of needs, is expressed as X ⁽⁰⁾; If image is RGB image three-colo(u)r, then image is transformed to YCbCr color space, Y-component is interpolated into the size of needs, is expressed as X ⁽⁰⁾if, X ⁽⁰⁾∈ R ^{n " × 1}, then the matrix of coefficients that A, B are N " × N " dimension is defined;

Second step: the value of design factor matrix A and B, it is divided into the following steps:

(1) by X ⁽⁰⁾being divided into size is image fritter, note become x _i, i=1,2 ... N, N represent the number of image block, have overlap between block adjacent when piecemeal; By X ⁽⁰⁾carry out convolution with a low frequency Gaussian convolution core, obtain low-frequency image, then by X ⁽⁰⁾subtract each other with low-frequency image, difference reflection X ⁽⁰⁾high-frequency information (here by its called after high frequency imaging), high frequency imaging being divided into size is image fritter, find and x _ithe image fritter of corresponding position, and comply with left-to-right, from top to bottom, form column vector by row reading manner, be designated as

(2) calculate with all { μ ₁, μ ₂... μ _kbetween Euclidean distance, find apart from minimum that, under it, be labeled as k _i, find { a ₁, a ₂... a _kin under be designated as k _iweighting parameter, as x _iautoregressive model weighting parameter

(3) center pixel value of establishing, χ _ifor x _imiddle x ' _ineighborhood territory pixel value composition vector, use formula (6) compute matrix A;

I, j are coordinate variables, get positive integer, and span is 1 ~ N ".

(4) at X ⁽⁰⁾in the fritter that image is divided into, find l similar fritter, variable l=1,2 ... L, L represent the quantity of similar fritter, are positive integer, represent X ⁽⁰⁾in and x _isimilar fritter, calculates in value, wherein represent normalized factor, h is constant, if for weight vector, represent the center pixel set of all similar fritters, then use formula (7) compute matrix B:

I, l are coordinate variables, get positive integer, and span is 1 ~ N ".

3rd step: set the constant γ in formula below ₁, γ ₂, γ ₃, P, e, Mid_Iter and maximum iterations Max_Iter, setting constant matrices τ, initialization iterations k=0;

4th step: establish I representation unit matrix, I matrix size and matrix A, B are the same, and D represents down-sampling matrix, D is according to the setting of reconstruction multiple, and H is Gaussian Blur matrix, and it is the matrix form of Gaussian convolution core, according to the setting of Gaussian convolution core, be a circular matrix, computing formula (8):

X ^(k+1/2)＝X ^(k)+γ ₃[(DH) ^TY-(DH) ^TDHX ^(k)]

-γ ₁(I-A) ^T(I-A)X ^(k)-γ ₂(I-B) ^T(I-B)X ^(k)(8)

X ^(k+1/2), X ^(k)represent the reconstructed results of iteration kth+1/2 and kth step respectively.

5th step: if R _irepresent x _ifrom X, cutting out, that is: x _i=R _ix, if iterations k is less than Mid_Iter, uses formula (9) compute sparse coefficient Θ ^(k+1/2), Θ ^(k+1/2)=[α ₁, α ₂... α _n]; Otherwise, use formula (10) to calculate α _i:

Θ ^(k+1/2)＝[ΨR ₁X ^(k+1/2)，ΨR ₂X ^(k+1/2)，…ΨR _NX ^(k+1/2)](9)

${\mathrm{\α}}_i=\underset{\mathrm{\α}}{\arg \min }\{{\left|\right|x_i-\mathrm{\Φ\α}\left|\right|}^2+{\mathrm{\γ}}_4{\left|\mathrm{\α}\right|}_1\}---\left(10\right)$

Wherein, γ ₄for constant, formula (10) adopts characteristic symbol finding algorithm to solve, and detailed process is as follows:

A () defines vectorial θ ∈ R ^{m × 1}, θ _ja jth element in representation vector θ, θ _j∈ {-1,0,1}, initialization definition of activities set β={ }, and be initialized as empty set;

B (), for the element in α being 0, calculates find out j value, α _jrepresent each element of jth of α, if then θ _j=-1, β=β ∪ j}, if: then θ _j=1, β=β ∪ { j};

C () selects the column vector composition being designated as β under in Φ select the element being designated as β under in α and θ to form respectively with calculate ${\stackrel{\‾}{\mathrm{\α}}}_{\mathrm{new}}={\left({\stackrel{\‾}{\mathrm{\Φ}}}^T\stackrel{\‾}{\mathrm{\Φ}}\right)}^{-1}({\stackrel{\‾}{\mathrm{\Φ}}}^T{x}_i-{\mathrm{\γ}}_4\stackrel{\‾}{\mathrm{\θ}}/2),$ Then check one by one with relevant position element, sees which element changes symbol (refer to by just become negative or by negative change just), is provided with the value reindexing of Num position, general the element of middle negate position sets to 0 (each Value Operations only to a position, the value of other position remains unchanged), uses the altered new vector of representative value, then there is Num kind value, will num kind value substitute into respectively in, find out and make be worth minimum that value, its assignment is given the element being designated as β under in α with middle relevant position element makes identical change, will the subscript becoming 0 with element in α removes from β, upgrades θ=sign (α);

(d) judge in α be not 0 element whether meet: $\frac{\∂{\left|\right|x_i-\mathrm{\Φ\α}\left|\right|}^2}{\∂{\mathrm{\α}}_j}+{\mathrm{\γ}}_4\mathrm{sign}\left({\mathrm{\α}}_j\right)=0,\∀{\mathrm{\α}}_j\≠0,$ If do not meet, then perform (c) step, otherwise judge in α be 0 element whether meet: if do not meet, perform (b) step, otherwise, return value (the i.e. α of α _i=α).

6th step: calculate Θ with formula (11) ^(k+1), defining operation (T _τ[Z]) _{i, j}=sign (Z _{i, j}) max{|Z _{i, j}|-τ _{i, j}, 0}, Z are the objects of this operation.

Θ ^(k+1)＝T _τ[Θ ^(k+1/2)](11)

7th step: use formula (12) to calculate X ^(k+1)

$X^{(k+1)}={\left(\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{R}_i^T{R}_i\right)}^{-1}\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{R}_i^T\mathrm{\Φ}{\mathrm{\α}}_i---\left(12\right)$

8th step: if mod (k, P)==0, and k>=Mid_Iter, by X ^(k+1)replace the X in second step ⁽⁰⁾, recalculate A and B, and use formula (13) to calculate τ _{i, j};

${\mathrm{\τ}}_{i,j}=\frac{c{\mathrm{\σ}}_n^2}{{\mathrm{\σ}}_{i,j}+\mathrm{\δ}}---\left(13\right)$

C is a constant, σ _nbe the standard deviation of picture noise, δ is a smaller constant, σ _{i, j}be calculated as follows: use X ^(k+1), extract fritter, and ask and x _isimilar fritter, to all and x _isimilar fritter vector calculate then propose l=1,2 ... the jth number of L, the standard deviation calculating these numbers is σ _{i, j}.

9th step: iterations k=k+1;

Tenth step: judge with k>=Max_Iter, wherein there is a condition to set up, then stop iteration, return super-resolution image X; Otherwise repeat the 4th step to the tenth step;

11 step: if input picture is gray level image, directly export X, if coloured image, is then interpolated into the size identical with X by CbCr component, then brightness X and color CbCr is converted into rgb space, exports the image after rebuilding.

The present invention is a kind of image super-resolution rebuilding method, compared with prior art, the invention has the advantages that image reconstruction better quality.

In order to the validity after verification algorithm reconstruction, by the same bi-cubic interpolation of result after reconstruction herein, [the JianchaoYang that JianchaoYang proposes, JohnWright, ThomasS.Huang, YiMa.Imagesuper-resolutionviasparserepresentation [J] .IEEETransactiononimageprocessing.2010, 19 (11): 2861-2873.] based on [ImageDeblurringandSuper-ResolutionbyAdaptiveSparseDomain SelectionandAdaptiveRegularization [J] .IEEETransactiononimageprocessing.2011 that method and the WeishengDong of rarefaction representation propose, 20 (7): 1838-1857.] ASDS reconstruction algorithm contrasts.

Accompanying drawing illustrates: Figure 1 shows that the visual effect after the reconstruction of four kinds of methods contrasts, and the image in Fig. 1 is through reducing process, and the upper left corner is the local original image intercepted.In Fig. 1 first is classified as the low-resolution image of input, second result being classified as bi-cubic interpolation, 3rd is classified as the result that JianchaoYang algorithm (be called for short Yang algorithm) rebuilds, 4th is classified as the result that WeishengDong algorithm (be called for short Dong algorithm) rebuilds, and the 5th is classified as the result that the inventive method rebuilds.Can see that the result that bi-cubic interpolation obtains is the fuzzyyest from result (the image upper left corner), and the Edge restoration degree of the image obtained based on the method for rarefaction representation is inadequate, image visual effect is also not as ASDS method, and ASDS rebuild after picture quality very high, but it is in butterfly's wing figure, cause partial distortion, method of the present invention can obtain best visual effect, does not almost have visible distortion.

In order to objective appraisal four kinds of super resolution ratio reconstruction methods, table one gives the data such as Y-PSNR (PSNR:PeakSignaltoNoiseRatio) and structural similarity index (SSIM:StructureSimilarityIndex) of four kinds of super resolution ratio reconstruction methods.Can obtain from table one, the value of PSNR and SSIM of bi-cubic interpolation is minimum, Yang algorithm is relative to bi-cubic interpolation, PSNR and SSIM value has and promotes largely, it is the highest that Dong algorithm and algorithm indices of the present invention improve degree, and algorithm of the present invention in most cases PSNR and SSIM index is slightly better than Dong algorithm.

Table one

Embodiment:

The technical solution used in the present invention is:

One, redundant dictionary, encoder dictionary parameter training:

for redundant dictionary, Ψ=[ψ ₁, ψ ₂..., ψ _n] ∈ R ^{m × n}for encoder dictionary, m, n are positive integer, wherein m=512, n=49, and double dictionary refers to redundant dictionary and encoder dictionary.

The first step: super-resolution image in reading images storehouse, transfers super-resolution image to gray level image, being then divided into size is fritter sample, the image fritter obtained is comply with from left to right, from top to bottom, by row reading manner formed column vector, use s _i∈ R ⁿ, i=1,2 ..., Q represents that the column vector that each fritter is formed, Q are the number of total column vector;

Second step: calculate s _ivariance Var (s _i), only retain Var (s _i) be greater than the vector of threshold value TH, wherein TH span is: 4.5 ~ 20, finally obtains training sample set S=[s ₁, s ₂... s _m], M is greater than 120000;

3rd step: formula (1) is solved, adopt alternative manner to solve redundant dictionary Φ, encoder dictionary Ψ, Θ are sparse coefficient, and λ, η are constant, and η gets the numerical value being approximately equal to 1, and λ span is 0.05 ~ 0.2, represent and ask l ₂norm, || ₁represent and ask l ₁norm:

$\mathrm{fuction}=\{\mathrm{\Φ},\mathrm{\Ψ},\mathrm{\Θ}\}=\arg \underset{\mathrm{\Φ},\mathrm{\Ψ},\mathrm{\Θ}}{\min }\{{\left|\right|S-\mathrm{\Φ\Θ}\left|\right|}_2^2+\mathrm{\η}{\left|\right|\mathrm{\Θ}-\mathrm{\ΨS}\left|\right|}_2^2+\mathrm{\λ}{\left|\mathrm{\Θ}\right|}_1\}---\left(1\right)$

(1) with gaussian random matrix initialisation redundant dictionary Φ, with unit matrix initialization codes dictionary Ψ, with full null matrix initialization sparse coefficient Θ, iterations k=0, maximum iteration time Max_Iter gets 800 ~ 1500, iteration convergence controlling elements ε=10 ^-6;

(2) (T is defined _ζ[O]) _{i, j}=sign (O _{i, j}) max{|O _{i, j}|-ζ, 0} are threshold operation operator, and ζ represents threshold operation variable, and O represents threshold operation matrix variables, O _{i, j}be designated as the element of (i, j) under in representing matrix O, sign () is symbol manipulation operator, gets σ _Θ=2|| Φ ^tΦ+η I|| _f, || || _frepresent and ask Frobenius norm, I is unit matrix, uses formula (2) to upgrade current Θ value;

${\mathrm{\Θ}}^{k+1}=T_{\mathrm{\λ}/2{\mathrm{\σ}}_{\mathrm{\Θ}}}[(1-\frac{\mathrm{\η}}{{\mathrm{\σ}}_{\mathrm{\Θ}}}){\mathrm{\Θ}}^k+\frac{1}{{\mathrm{\σ}}_{\mathrm{\Θ}}}({\mathrm{\Φ}}^T(S-\mathrm{\Φ}{\mathrm{\Θ}}^k)+\mathrm{\η\ΨS})]---\left(2\right)$

Wherein Θ ^k+1, Θ ^krepresent the Θ value of iteration kth+1 and kth step respectively, Φ ^trepresent the transposition of Φ.

(3) defining operation π (d)=d/max (1, || d||), d is vector, and this operation represents by vector projection to unit length, definition σ _Φ=2|| Θ Θ ^t|| _f, use formula (3) to upgrade current Φ value:

${\mathrm{\Φ}}^{k+1}=\mathrm{\π}({\mathrm{\Φ}}^k+\frac{1}{{\mathrm{\σ}}_{\mathrm{\Φ}}}(S-{\mathrm{\Φ}}^k\mathrm{\Θ}){\mathrm{\Θ}}^T)---\left(3\right)$

π () expression herein carries out unit length projection, wherein Φ to each row of Φ ^k+1, Φ ^krepresent the Φ value of iteration kth+1 and kth step respectively, Θ ^trepresent the transposition of Θ.

(4) σ is calculated _Ψ=2||SS ^t|| _f, use formula (4) to upgrade current Ψ value:

${\mathrm{\Ψ}}^{k+1}=\mathrm{\π}({\mathrm{\Ψ}}^k+\frac{1}{{\mathrm{\σ}}_{\mathrm{\Ψ}}}(\mathrm{\Θ}-{\mathrm{\Ψ}}^{k}S)S^T)---\left(4\right)$

π () expression herein carries out unit length projection to each row of Ψ, wherein Ψ ^k+1, Ψ ^krepresent the Ψ value of iteration kth+1 and kth step respectively, S ^trepresent the transposition of S.

(5) iterations k=k+1;

(6) Θ, Ψ, Φ value of current calculating and the front value once calculated are substituted into formula (1), the value of calculating target function respectively, judge ${\left|\right|{\mathrm{function}}^{(k+1)}-{\mathrm{function}}^{\left(k\right)}\left|\right|}_2^2/{\left|\right|{\mathrm{function}}^{\left(k\right)}\left|\right|}_2^2<\mathrm{\&epsiv;}$ (function ^(k+1), function ^krepresent the function value that kth+1 and kth step calculate respectively) and k>=Max_Iter condition whether meet, meet wherein any one, stopping iteration, output Φ value and Ψ value; Otherwise repeat (2) to (5).

Two, autoregressive model weighting parameter training:

The first step: super-resolution image in reading images storehouse, transfer image to gray level image, then with one low frequency Gaussian convolution core carries out convolution, and (Gaussian convolution core size is 7 × 7, standard deviation is 1.6), obtain low-frequency image, then original image and low-frequency image are subtracted each other, the high-frequency information (here by its called after high frequency imaging) of difference reflection image, high frequency imaging being divided into size is image fritter, n=49, finds and s in super-resolution image _ithe image fritter of corresponding position, and comply with left-to-right, from top to bottom, form column vector by row reading manner, be designated as all S=[s ₁, s ₂... s _m] corresponding high frequency fritter vector set is combined into $S^h=[s_1^h,s_2^h,...s_M^h].$

Second step: use K-means sorting algorithm by S ^hbe divided into K class { C ₁, C ₂... C _k, k gets 200, m _krepresent C _kthe number of middle vector, uses formula (5) to calculate the barycenter μ of each class _k, k=1,2 ..., K, according to S ^hclassification results, by S=[s ₁, s ₂... s _m] be also divided into K class, be expressed as { S ₁, S ₂... S _k}:

$\begin{array}{cc}{\mathrm{\&mu;}}_k=\frac{1}{m_k}\underset{i=1}{\overset{m_k}{\mathrm{\&Sigma;}}}{s}_i^h& s_i^h\&Element;C_k\end{array}---\left(5\right)$

3rd step: if s ' _irepresent s _icenter pixel value, q _ifor s _imiddle s ' _ineighborhood territory pixel value composition vector, Size of Neighborhood is got 3 × 3 and (is comprised center pixel, q _ifor the vector formed after removing center pixel, be the column vector of 8 elements), use least square method calculating formula in a _k, a _kbe the vector of 8 × 1, identical process is carried out to all classification, obtain autoregressive model weighting parameter combination { a ₁, a ₂... a _k;

4th step: output from regression model weighting parameter combination { a ₁, a ₂... a _kand the barycenter { μ of each class ₁, μ ₂... μ _k;

Three, image super-resolution rebuilding

Redundant dictionary Φ, encoder dictionary Ψ, autoregressive model weighting parameter combination { a ₁, a ₂... a _k, the barycenter { μ of each class ₁, μ ₂... μ _kfor training in advance obtains, once training can use always.

The first step: read in the low-resolution image Y needing to carry out rebuilding, if gray level image, uses bi-cubic interpolation Y to be interpolated into the size of needs, is expressed as X ⁽⁰⁾; If image is RGB image three-colo(u)r, then image is transformed to YCbCr color space, Y-component is interpolated into the size of needs, is expressed as X ⁽⁰⁾if, X ⁽⁰⁾∈ R ^{n " × 1}, then the matrix of coefficients that A, B are N " × N " dimension is defined; ;

Second step: design factor matrix A and B, it is divided into the following steps:

(1) by X ⁽⁰⁾being divided into size is image fritter, note become x _i, i=1,2 ... N, N represent the number (size variation with input picture) of image block, have overlap (transverse direction or longitudinal overlap 4 pixel width) between block adjacent when piecemeal; By X ⁽⁰⁾carry out convolution with a low frequency Gaussian convolution core, obtain low-frequency image, then by X ⁽⁰⁾subtract each other with low-frequency image, difference reflection X ⁽⁰⁾high-frequency information (here by its called after high frequency imaging), high frequency imaging being divided into size is image fritter, find and x _ithe image fritter of corresponding position, and comply with left-to-right, from top to bottom, form column vector by row reading manner, be designated as

(2) calculate with all { μ ₁, μ ₂... μ _kbetween Euclidean distance, find apart from minimum that, under it, be labeled as k _i, find { a ₁, a ₂... a _kin under be designated as k _iweighting parameter, as x _iautoregressive model weighting parameter

(3) if x ' _irepresent x _icenter pixel value, χ _ifor x _imiddle x ' _ineighborhood territory pixel value composition vector, use formula (6) compute matrix A;

I, j are coordinate variables, get positive integer, and span is 1 ~ N ".

(4) at X ⁽⁰⁾in the fritter be divided in image, find l similar fritter, variable l=1,2 ... L, L represent the quantity of similar fritter, are positive integer, and L span is 7 ~ 10, represent X ⁽⁰⁾in and x _isimilar fritter, calculates in value, wherein represent normalized factor, h is constant, and span is 65 ~ 70, if for weight vector, represent the center pixel set of all similar fritters, then use formula (7) compute matrix B:

I, l are coordinate variables, get positive integer, and span is 1 ~ N ".

3rd step: presetting: γ ₁span 0.008 ~ 0.01, γ ₂span 0.04 ~ 0.1, γ ₃get the value of about 6.5, P=20, e=10 ^-6, Mid_Iter=100 and maximum iterations Max_Iter=150, setting constant matrices τ=0, initialization iterations k=0;

4th step: establish I representation unit matrix, I matrix size and matrix A, B is the same, D represents down-sampling matrix, D is according to the setting of reconstruction multiple, H is Gaussian Blur matrix, it is that Gaussian convolution core is (when the reconstruction factor is 3, Gaussian convolution core size is 7 × 7, standard deviation is 1.6, when the reconstruction factor is 2, Gaussian convolution core size is 5 × 5, standard deviation is about 0.9 ~ 1.1, when the reconstruction factor is 4, Gaussian convolution core size is 7 × 7, standard deviation is 1.7 ~ 1.8) matrix form, set according to Gaussian convolution core, it is a circular matrix, computing formula (8):

X ^(k+1/2)＝X ^(k)+γ ₃[(DH) ^TY-(DH) ^TDHX ^(k)]

-γ ₁(I-A) ^T(I-A)X ^(k)-γ ₂(I-B) ^T(I-B)X ^(k)(8)

X ^(k+1/2), X ^(k)represent the reconstructed results of iteration kth+1/2 and kth step respectively.

5th step: if R _irepresent x _ifrom X, cutting out, that is: x _i=R _ix, if iterations k is less than Mid_Iter, uses formula (9) compute sparse coefficient Θ ^(k+1/2), Θ ^(k+1/2)=[α ₁, α ₂... α _n]; Otherwise, use formula (10) to calculate α _i:

Θ ^(k+1/2)＝[ΨR ₁X ^(k+1/2)，ΨR ₂X ^(k+1/2)，…ΨR _NX ^(k+1/2)](9)

${\mathrm{\&alpha;}}_i=\underset{\mathrm{\&alpha;}}{\arg \min }\{{\left|\right|x_i-\mathrm{\&Phi;\&alpha;}\left|\right|}^2+{\mathrm{\&gamma;}}_4{\left|\mathrm{\&alpha;}\right|}_1\}---\left(10\right)$

Wherein, γ ₄for constant, span is 0.1 ~ 0.2, and formula (10) adopts characteristic symbol finding algorithm to solve, and detailed process is as follows:

A () defines vectorial θ ∈ R ^{m × 1}, θ _ja jth element in representation vector θ, θ _j∈ {-1,0,1}, initialization definition of activities set β={ }, and be initialized as empty set;

B (), for the element in α being 0, calculates find out j value, α _jrepresent each element of jth of α, if then θ _j=-1, β=β ∪ j}, if: then θ _j=1, β=β ∪ { j};

C () selects the column vector composition being designated as β under in Φ select the element being designated as β under in α and θ to form respectively with calculate ${\stackrel{\&OverBar;}{\mathrm{\&alpha;}}}_{\mathrm{new}}={\left({\stackrel{\&OverBar;}{\mathrm{\&Phi;}}}^T\stackrel{\&OverBar;}{\mathrm{\&Phi;}}\right)}^{-1}({\stackrel{\&OverBar;}{\mathrm{\&Phi;}}}^T{x}_i-{\mathrm{\&gamma;}}_4\stackrel{\&OverBar;}{\mathrm{\&theta;}}/2),$ Then check one by one with relevant position element, sees which element changes symbol (refer to by just become negative or by negative change just), is provided with the value reindexing of Num position, general the element of middle negate position sets to 0 (each Value Operations only to a position, the value of other position remains unchanged), uses the altered new vector of representative value, then there is Num kind value, will num kind value substitute into respectively in, find out and make be worth minimum that value, its assignment is given the element being designated as β under in α with middle relevant position element makes identical change, will the subscript becoming 0 with element in α removes from β, upgrades θ=sign (α);

(d) judge in α be not 0 element whether meet: $\frac{\&PartialD;{\left|\right|x_i-\mathrm{\&Phi;\&alpha;}\left|\right|}^2}{\&PartialD;{\mathrm{\&alpha;}}_j}+{\mathrm{\&gamma;}}_4\mathrm{sign}\left({\mathrm{\&alpha;}}_j\right)=0,\&ForAll;{\mathrm{\&alpha;}}_j\&NotEqual;0,$ If do not meet, then perform (c) step, otherwise judge in α be 0 element whether meet: if do not meet, perform (b) step, otherwise, return value (the i.e. α of α _i=α);

6th step: calculate Θ with formula (11) ^(k+1), defining operation (T _τ[Z]) _{i, j}=sign (Z _{i, j}) max{|Z _{i, j}|-τ _{i, j}, 0}, Z are the objects of this operation:

Θ ^(k+1)＝T _τ[Θ ^(k+1/2)](11)

7th step: use formula (12) to calculate X ^(k+1):

$X^{(k+1)}={\left(\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{R}_i^T{R}_i\right)}^{-1}\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{R}_i^T\mathrm{\&Phi;}{\mathrm{\&alpha;}}_i---\left(12\right)$

8th step: if mod (k, P)==0, and k>=Mid_Iter, by X ^(k+1)replace the X in second step ⁽⁰⁾, recalculate A and B, and use formula (13) to calculate τ _{i, j}:

${\mathrm{\&tau;}}_{i,j}=\frac{c{\mathrm{\&sigma;}}_n^2}{{\mathrm{\&sigma;}}_{i,j}+\mathrm{\&delta;}}---\left(13\right)$

C is a constant, σ _nthe standard deviation of picture noise, c σ _nspan is 0.1 ~ 3.6, and normal image gets 0.1 ~ 0.6; δ is a smaller constant, gets 0.35, σ _{i, j}be calculated as follows: use X ^(k+1), extract fritter, and ask and x _isimilar fritter, to all and x _isimilar fritter vector calculate then propose l=1,2 ... the jth number of L, the standard deviation calculating these numbers is σ _{i, j}.

9th step: iterations k=k+1;

Tenth step: judge with k>=Max_Iter, wherein there is a condition to set up, then stop iteration, return super-resolution image X; Otherwise repeat the 4th step to the tenth step;

11 step: if input picture is gray level image, directly export X, if coloured image, is then interpolated into the size identical with X by CbCr component, then brightness X and color CbCr is converted into rgb space, exports the image after rebuilding.

Claims (1)

1., based on the compressed sensing image super-resolution rebuilding method of double dictionary study, it is characterized in that following steps:

A, redundant dictionary, encoder dictionary parameter training:

for redundant dictionary, Ψ=[ψ ₁, ψ ₂..., ψ _n] ∈ R ^{m × n}for encoder dictionary, m, n are positive integer, wherein m=512, n=49, and double dictionary refers to and produces redundant dictionary and encoder dictionary two dictionaries simultaneously;

The first step: super-resolution image in reading images storehouse, transfers super-resolution image to gray level image, being then divided into size is fritter sample, by the image fritter that obtains according to from left to right, from top to bottom, form column vector by row reading manner, use s _i∈ R ⁿ, i=1,2 ..., Q represents that the column vector that each fritter is formed, Q are the number of total column vector;

Second step: calculate s _ivariance Var (s _i), only retain Var (s _i) be greater than the vector of threshold value TH, wherein TH span is: 4.5 ~ 20, finally obtains training sample set S=[s ₁, s ₂... s _m], M is greater than 120000;

3rd step: formula (1) is solved, adopt alternative manner to solve redundant dictionary Φ, encoder dictionary Ψ, Θ are sparse coefficient, and λ, η are constant, and η gets the numerical value equaling 1, and λ span is 0.05 ~ 0.2, represent and ask l ₂norm, || ₁represent and ask l ₁norm:

$\mathrm{fuction}=\{\mathrm{\&Phi;},\mathrm{\&Psi;},\mathrm{\&Theta;}\}=\arg \underset{\mathrm{\&Phi;},\mathrm{\&Psi;},\mathrm{\&Theta;}}{\min }\{{\left|\right|S-\mathrm{\&Phi;\&Theta;}\left|\right|}_2^2+\mathrm{\&eta;}{\left|\right|\mathrm{\&Theta;}-\mathrm{\&Psi;S}\left|\right|}_2^2+\mathrm{\&lambda;}{\left|\mathrm{\&Theta;}\right|}_1\}---\left(1\right)$

(1) with gaussian random matrix initialisation redundant dictionary Φ, with unit matrix initialization codes dictionary Ψ, with full null matrix initialization sparse coefficient Θ, iterations k=0, maximum iteration time Max_Iter gets 800 ~ 1500, iteration convergence controlling elements ε=10 ^-6;

(2) (T is defined _ζ[O]) _{i, j}=sign (O _{i, j}) max{|O _{i, j}|-ζ, 0} are threshold operation operator, and ζ represents threshold operation variable, and O represents threshold operation matrix variables, O _{i, j}be designated as the element of (i, j) under in representing matrix O, sign () is symbol manipulation operator, gets σ _Θ=2|| Φ ^tΦ+η I|| _f, || || _frepresent and ask Frobenius norm, I is unit matrix, uses formula (2) to upgrade current Θ value;

${\mathrm{\&Theta;}}^{k+1}=T_{\mathrm{\&lambda;}/2{\mathrm{\&sigma;}}_{\mathrm{\&Theta;}}}[(1-\frac{\mathrm{\&eta;}}{{\mathrm{\&sigma;}}_{\mathrm{\&Theta;}}}){\mathrm{\&Theta;}}^k+\frac{1}{{\mathrm{\&sigma;}}_{\mathrm{\&Theta;}}}({\mathrm{\&Phi;}}^T(S-\mathrm{\&Phi;}{\mathrm{\&Theta;}}^k)+\mathrm{\&eta;\&Psi;S})]---\left(2\right)$

Wherein Θ ^k+1, Θ ^krepresent the Θ value of iteration kth+1 and kth step respectively, Φ ^trepresent the transposition of Φ;

(3) defining operation π (d)=d/max (1, || d||), d is vector, and this operation represents by vector projection to unit length, definition σ _Φ=2|| Θ Θ ^t|| _f, use formula (3) to upgrade current Φ value:

${\mathrm{\&Phi;}}^{k+1}=\mathrm{\&pi;}({\mathrm{\&Phi;}}^k+\frac{1}{{\mathrm{\&sigma;}}_{\mathrm{\&Phi;}}}(S-{\mathrm{\&Phi;}}^k\mathrm{\&Theta;}){\mathrm{\&Theta;}}^T)---\left(3\right)$

π () expression herein carries out unit length projection, wherein Φ to each row of Φ ^k+1, Φ ^krepresent the Φ value of iteration kth+1 and kth step respectively, Θ ^trepresent the transposition of Θ;

(4) σ is calculated _Ψ=2||SS ^t|| _f, use formula (4) to upgrade current Ψ value:

${\mathrm{\&Psi;}}^{k+1}=\mathrm{\&pi;}({\mathrm{\&Psi;}}^k+\frac{1}{{\mathrm{\&sigma;}}_{\mathrm{\&Psi;}}}(\mathrm{\&Theta;}-{\mathrm{\&Psi;}}^{k}S)S^T)---\left(4\right)$

π () expression herein carries out unit length projection to each row of Ψ, wherein Ψ ^k+1, Ψ ^krepresent the Ψ value of iteration kth+1 and kth step respectively, S ^trepresent the transposition of S;

(5) iterations k=k+1;

(6) Θ, Ψ, Φ value of current calculating and the front value once calculated are substituted into formula (1), the value of calculating target function respectively, judge whether meet with k>=Max_Iter condition, meet wherein any one, stop iteration, export Φ value and Ψ value, otherwise repeat (2) to (5), wherein function ^(k+1), function ^krepresent the function value of kth+1 and the calculating of kth step respectively;

B, autoregressive model weighting parameter are trained:

The first step: super-resolution image in reading images storehouse, transfer image to gray level image, then with one size is 7 × 7, standard deviation be 1.6 low frequency Gaussian convolution core carry out convolution, obtain low-frequency image, then original image and low-frequency image are subtracted each other, the high-frequency information of difference reflection image, we are by its called after high frequency imaging, high frequency imaging are divided into size and are image fritter, n=49, finds and s in super-resolution image _ithe image fritter of corresponding position, and according to from left to right, from top to bottom, form column vector by row reading manner, be designated as all S=[s ₁, s ₂... s _m] corresponding high frequency fritter vector set is combined into

Second step: use K-means sorting algorithm by S ^hbe divided into K class { C ₁, C ₂... C _k, k gets 200, m _krepresent C _kthe number of middle vector, uses formula (5) to calculate the barycenter μ of each class _k, k=1,2 ..., K, according to S ^hclassification results, by S=[s ₁, s ₂... s _m] be also divided into K class, be expressed as { S ₁, S ₂... S _k}:

${\mathrm{\&mu;}}_k=\frac{1}{m_k}\underset{i=1}{\overset{m_k}{\mathrm{\&Sigma;}}}{s}_i^h,s_i^h\&Element;C_k---\left(5\right)$

3rd step: if s ' _irepresent s _icenter pixel value, q _ifor s _imiddle s ' _ineighborhood territory pixel value composition vector, Size of Neighborhood gets 3 × 3, q _ifor the vector formed after removing center pixel, be the column vector of 8 elements, use least square method calculating formula in a _k, a _kbe the vector of 8 × 1, identical process is carried out to all classification, obtain autoregressive model weighting parameter combination { a ₁, a ₂... a _k;

4th step: output from regression model weighting parameter combination { a ₁, a ₂... a _kand the barycenter { μ of each class ₁, μ ₂... μ _k;

C, image super-resolution rebuilding

Redundant dictionary Φ, encoder dictionary Ψ, autoregressive model weighting parameter combination { a ₁, a ₂... a _k, the barycenter { μ of each class ₁, μ ₂... μ _kfor training in advance obtains, once training can use always;

The first step: read in the low-resolution image Y needing to carry out rebuilding, if gray level image, uses bi-cubic interpolation Y to be interpolated into the size of needs, is expressed as X ⁽⁰⁾; If image is RGB image three-colo(u)r, then image is transformed to YCbCr color space, Y-component is interpolated into the size of needs, is expressed as X ⁽⁰⁾if, then define the matrix of coefficients that A, B are N " × N " dimension;

Second step: design factor matrix A and B, it is divided into the following steps:

(1) by X ⁽⁰⁾being divided into size is image fritter, note become x _i, i=1,2 ... N, N represent the number of image block, and the value of N, with the size variation of input picture, has overlap between block adjacent when piecemeal, transverse direction or longitudinal overlap 4 pixel width; By X ⁽⁰⁾carry out convolution with a low frequency Gaussian convolution core, obtain low-frequency image, then by X ⁽⁰⁾subtract each other with low-frequency image, difference reflection X ⁽⁰⁾high-frequency information, here by its called after high frequency imaging, high frequency imaging being divided into size is image fritter, find and x _ithe image fritter of corresponding position, and according to from left to right, from top to bottom, form column vector by row reading manner, be designated as

(2) calculate with all { μ ₁, μ ₂... μ _kbetween Euclidean distance, find apart from minimum that, under it, be labeled as k _i, find { a ₁, a ₂... a _kin under be designated as k _iweighting parameter, as x _iautoregressive model weighting parameter

(3) if x ' _irepresent x _icenter pixel value, χ _ifor x _imiddle x ' _ineighborhood territory pixel value composition vector, use formula (6) compute matrix A;

I, j are coordinate variables, get positive integer, and span is 1 ~ N ";

(4) at X ⁽⁰⁾in the fritter be divided in image, find l similar fritter, variable l=1,2 ... L, L represent the quantity of similar fritter, are positive integer, and L span is 7 ~ 10, represent X ⁽⁰⁾in and x _isimilar fritter, calculates in value, wherein represent normalized factor, h is constant, and span is 65 ~ 70, if for weight vector, represent the center pixel set of all similar fritters, then use formula (7) compute matrix B:

I, l are coordinate variables, get positive integer, and span is 1 ~ N ";

3rd step: presetting: γ ₁span 0.008 ~ 0.01, γ ₂span 0.04 ~ 0.1, γ ₃value 6.5, P=20, e=10 ^-6, Mid_Iter=100 and maximum iterations Max_Iter=150, setting constant matrices τ=0, initialization iterations k=0;

4th step: establish I representation unit matrix, I matrix size and matrix A, B are the same, D represents down-sampling matrix, D is according to the setting of reconstruction multiple, and H is Gaussian Blur matrix, and it is Gaussian convolution core, when the reconstruction factor is 3, Gaussian convolution core size is 7 × 7, and standard deviation is 1.6, when the reconstruction factor is 2, Gaussian convolution core size is 5 × 5, standard deviation is 0.9 ~ 1.1, and when the reconstruction factor is 4, Gaussian convolution core size is 7 × 7, standard deviation is 1.7 ~ 1.8, according to the setting of Gaussian convolution core, reconstruction matrix is a circular matrix, computing formula (8):

X ^(k+1/2)＝X ^(k)+γ ₃[(DH) ^TY-(DH) ^TDHX ^(k)]

-γ ₁(I-A) ^T(I-A)X ^(k)-γ ₂(I-B) ^T(I-B)X ^(k)(8)

X ^(k+1/2), X ^(k)represent the reconstructed results of iteration kth+1/2 and kth step respectively;

5th step: if R _irepresent x _ifrom X, cutting out, that is: x _i=R _ix, if iterations k is less than Mid_Iter, uses formula (9) compute sparse coefficient Θ ^(k+1/2), Θ ^(k+1/2)=[α ₁, α ₂... α _n]; Otherwise, use formula (10) to calculate α _i:

Θ ^(k+1/2)＝[ΨR ₁X ^(k+1/2)，ΨR ₂X ^(k+1/2)，…ΨR _NX ^(k+1/2)](9)

${\mathrm{\&alpha;}}_i=\underset{\mathrm{\&alpha;}}{\arg \min }\{{\left|\right|x_i-\mathrm{\&Phi;\&alpha;}\left|\right|}^2+{\mathrm{\&gamma;}}_4{\left|\mathrm{\&alpha;}\right|}_1\}---\left(10\right)$

Wherein, γ ₄for constant, span is 0.1 ~ 0.2, and formula (10) adopts characteristic symbol finding algorithm to solve, and detailed process is as follows:

A () defines vectorial θ ∈ R ^{m × 1}, θ _ja jth element in representation vector θ, θ _j∈ {-1,0,1}, initialization definition of activities set β={ }, and be initialized as empty set;

B (), for the element in α being 0, calculates find out j value, α _jrepresent each element of jth of α, if then θ _j=-1, β=β ∪ j}, if: then θ _j=1, β=β ∪ { j};

C () selects the column vector composition being designated as β under in Φ select the element being designated as β under in α and θ to form respectively with calculate then check one by one with relevant position element, sees which element changes symbol, and sign modification refers to by just becoming negative or by negative change just, being provided with the value reindexing of Num position, will the element of middle negate position sets to 0, each Value Operations only to a position, and the value of other position remains unchanged, and uses the altered new vector of representative value, then there is Num kind value, will num kind value substitute into respectively in, find out and make be worth minimum that value, its assignment is given the element being designated as β under in α with middle relevant position element makes identical change, will the subscript becoming 0 with element in α removes from β, upgrades θ=sign (α);

(d) judge in α be not 0 element whether meet: if do not meet, then perform (c) step, otherwise judge in α be 0 element whether meet: if do not meet, perform (b) step, otherwise, return the value of α;

6th step: calculate Θ with formula (11) ^(k+1), defining operation (T _τ[Z]) _{i, j}=sign (Z _{i, j}) max{|Z _{i, j}|-τ _{i, j}, 0}, Z are the objects of this operation:

Θ ^(k+1)＝T _τ[Θ ^(k+1/2)](11)

7th step: use formula (12) to calculate X ^(k+1):

$X^{(k+1)}={\left(\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{R}_i^T{R}_i\right)}^{-1}\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{R}_i^T\mathrm{\&Phi;}{\mathrm{\&alpha;}}_i---\left(12\right)$

8th step: if mod (k, P)==0, and k>=Mid_Iter, by X ^(k+1)replace the X in second step ⁽⁰⁾, recalculate A and B, and use formula (13) to calculate τ _{i, j}:

${\mathrm{\&tau;}}_{i,j}=\frac{c{\mathrm{\&sigma;}}_n^2}{{\mathrm{\&sigma;}}_{i,j}+\mathrm{\&delta;}}---\left(13\right)$

C is a constant, σ _nthe standard deviation of picture noise, c σ _nspan is 0.1 ~ 3.6, and normal image gets 0.1 ~ 0.6; δ is a smaller constant, gets 0.35, σ _{i, j}be calculated as follows: use X ^(k+1), extract fritter, and ask and x _isimilar fritter, to all and x _isimilar fritter vector calculate then propose jth number, calculate these number standard deviations be σ _{i, j};

9th step: iterations k=k+1;

Tenth step: judge with k>=Max_Iter, wherein there is a condition to set up, then stop iteration, return super-resolution image X; Otherwise repeat the 4th step to the tenth step;

11 step: if input picture is gray level image, directly export X, if coloured image, is then interpolated into the size identical with X by CbCr component, then brightness X and color CbCr is converted into rgb space, exports the image after rebuilding.