CN104700436A - Edge constraint based image reconstruction method under multivariate observation - Google Patents

Abstract

The invention discloses an edge constraint based image reconstruction method under multivariate observation and mainly aims to solve the problems of the prior art of compressed sensing image reconstruction inaccuracy and low robustness. The edge constraint based image reconstruction method includes: 1) receiving an observation matrix, a low frequency wavelet decomposition coefficient and a multivariate measurement matrix; 2) acquiring a nonzero coefficient group supporting set through edge detection and relevant guides; 3) reconstructing high frequency wavelet coefficient in the nonzero coefficient group supporting set on the basis of a multivariate Gaussian model according to the observation matrix, the multivariate measurement matrix, basic covariance and residual covariance matrix in the Gibbs sampling method; 4) converting the low frequency wavelet decomposition coefficient and the reconstructed high frequency wavelet coefficient to obtain reconstruction images. Compared with OMP and BEPA method, the edge constraint based image reconstruction method has the advantages of high reconstruction image quality and good robustness, and can be reconstruction of natural images and medical images.

Description

Based on the image reconstructing method of edge constraint under multivariate observation

Technical field

The invention belongs to technical field of image processing, be specifically related to statistics compressed sensing image reconstructing method, can be used for being reconstructed natural image.

Background technology

In recent years, a kind of new data theory compressed sensing CS has been there is in signal transacting field, this theory realizes compression while data acquisition, breach tradition how Kui gather the restriction of this special sampling thheorem, for data acquisition technology brings revolutionary change, this theory is had broad application prospects in fields such as compression imaging system, military cryptology, wireless sensings.Compressive sensing theory mainly comprises three aspects such as reconstruct of the rarefaction representation of signal, the observation of signal and signal.Wherein designing restructing algorithm is fast and effectively CS theory successfully promoted and be applied to the important step of real data model and acquisition system.

In compression sampling field, wavelet basis is one group of good sparse base.The coefficient of dissociation that image obtains after wavelet decomposition, be divided into low frequency part and HFS, the low frequency that low frequency part comprises original image is sparse, it has been generally acknowledged that right and wrong are sparse, and HFS comprise image level, vertical, to angle information, have good openness.At present, often adopt and under wavelet field, low frequency is all retained, the method for sampling compressing observation is carried out to high frequency.The advantage of this method of sampling effectively can improve reconstructed image quality.

The people such as Lihan He propose the Bayes's compressed sensing image reconstructing method based on wavelet tree structure in document " Exploiting Structure in Wavelet-Based Bayesian Compressive Sensing ".The method sets up single Gauss model to multi-scale wavelet coefficient hierarchical, and by gibbs sampler reconstructed image.But image spread becomes column vector to carry out observation reconstruct by the method, not only not in conjunction with the priori of raw image data, and requires very high to calculator memory, limits the size of process image.

The people such as Jiao Wu propose the multivariable compression perceptual image reconstruct MPA based on mixed-scale model in document " Multivariate Compressive Sensing for Image Reconstruction in the Wavelet Domain:Using Scale Mixture Models ".The method is to wavelet coefficient structure multivariate distributed model, wavelet coefficient is caught to have this feature of aggregation, modeling is carried out to its statistic correlation, but the method have ignored the wavelet low frequency coefficient that remains to the directive function of Image Reconstruction, thus cause it not have robustness, and the image reconstructed is not accurate enough.

Summary of the invention

The object of the invention is to the deficiency of the above-mentioned prior art of pin, propose a kind of under multivariate observation based on the image reconstructing method of edge constraint, to make full use of the directive function of low-frequency wavelet coefficients to Image Reconstruction of reservation, improve the accuracy of reconstructed image.

Existing the object of the invention technical thought is: by setting up multivariate Gaussian model to multivariate calculation matrix, catch the aggregation of small echo; By the associating of rim detection and correlativity, the non-zero determining wavelet coefficient is instructed to support; By utilizing Gibbs sampling method iteration renewal successively to non-zero support coefficient, realize high-quality compressed sensing Image Reconstruction.

According to above-mentioned thinking, technical scheme of the present invention comprises the steps:

1. under multivariate observation based on the image reconstructing method of edge constraint, comprise the steps:

(1) take over party receives Random Orthogonal Gauss observing matrix Φ, low frequency wavelet coefficient of dissociation L, the horizontal high-frequent subband multivariate calculation matrix Y of image transmit leg transmission ₁, vertical high frequency subband multivariate calculation matrix Y ₂with diagonal angle high-frequency sub-band multivariate calculation matrix Y ₃, three high-frequency sub-band multivariate calculation matrix unification Y are represented;

(2) according to the observing matrix Φ, the low frequency wavelet coefficient of dissociation L that receive and high-frequency sub-band multivariate calculation matrix Y, the set of nonzero coefficient group index is obtained by the guidance of rim detection and correlativity: u={s ₁, s ₂..., s _i..., s _c, wherein s _irepresent the index of i-th nonzero coefficient group, i=1,2 ..., c, c are the columns being less than Φ:

(2.1) high-frequency sub-band being all zero by the low frequency wavelet coefficient of dissociation L and three of reception carries out wavelet inverse transformation, obtains Edge obscuring image;

(2.2) edge blurred picture carries out rim detection, obtains marginal position;

(2.3) pixel extracting position, corresponding edge in Edge obscuring image obtains fuzzy edge;

(2.4) one deck wavelet transformation is carried out to fuzzy edge, obtain fuzzy edge small echo high frequency coefficient and fuzzy edge wavelet low frequency coefficient; The position that fuzzy edge small echo high frequency coefficient absolute value is greater than threshold value h is set to 1, the position being less than threshold value h is set to 0 and obtains initial fuzzy location matrix, initial fuzzy location matrix is arranged in the Multivariable Fuzzy location matrix E of M × Q dimension according to the form of multi-variable matrix, wherein M is the columns of Φ, Q is the columns of Y, threshold value h=0.2;

(2.5) to be multiplied with high-frequency sub-band multivariate calculation matrix Y the correlation matrix Φ obtained according to the transposition of observing matrix Φ ^t* Y, by the absolute value of this correlation matrix | Φ ^t* Y| and Multivariable Fuzzy location matrix E is weighted summation, obtain coefficient importance matrix V=| Φ ^t* Y|+w*E, wherein w is weighting coefficient;

(2.6) every a line of coefficient importance matrix V is added, obtains coefficient importance vector, the index of c maximum in a coefficient importance vector element is taken out and forms nonzero coefficient group index set u={s ₁, s ₂..., s _i..., s _c, wherein i=1,2 ..., c, s _irepresent the index of i-th nonzero coefficient group, c is the integer being less than M;

(3) initialization outer iteration number of times t=1, model solution iteration total degree is N ₁, it is N that coefficient solves iteration total degree ₂; Reconstruct small echo high frequency coefficient matrix is treated in initialization be the null matrix of M × Q dimension with total stack result O, wherein M is the columns of Φ, and Q is the columns of Y; Reconstruct small echo high frequency coefficient matrix is treated in initialization every a line all obey basic covariance matrix Ω multivariate Gaussian distribution; Initialization residual error every a line all obey residual covariance matrix Π multivariate Gaussian distribution;

(4) according to observing matrix Φ, high-frequency sub-band multivariate calculation matrix Y, basic covariance matrix Ω, residual covariance matrix Π and nonzero coefficient group index set u, calculated by Gibbs sampling method and treat reconstruct small echo high frequency coefficient matrix

(4.1) initialization iterations i=1;

(4.2) obtain treating reconstruct small echo high frequency coefficient matrix according to basic covariance matrix Ω, residual covariance matrix Π, observing matrix Φ s _ithe covariance matrix of row:

Wherein for the s of observing matrix Φ _irow, for transposition;

(4.3) according to residual covariance matrix Π, observing matrix Φ with treat reconstruct small echo high frequency coefficient matrix s _ithe covariance matrix of row obtain treating reconstruct small echo high frequency coefficient matrix s _ithe mean vector of row:

Wherein for the kth of observing matrix Φ arranges, x _kfor treating reconstruct small echo high frequency coefficient matrix row k, wherein k is less than or equal to M for being more than or equal to 1 and being not equal to s _iinteger;

(4.4) according to treating reconstruct small echo high frequency coefficient matrix s _ithe mean vector of row and covariance matrix set up corresponding multivariate Gaussian model, calculate and treat reconstruct small echo high frequency coefficient matrix s _irow coefficient $x_{s_i}=\mathrm{Gaussian}({\mathrm{\μ}}_{s_i},{\mathrm{\Λ}}_{s_i}),$ Wherein, represent that stochastic generation one is obeyed mean vector and is covariance matrix is multivariate Gaussian distribution vector;

(4.5) iterations i and nonzero coefficient group index set u element number c is compared: if iterations i<c, then iterations i is from adding 1, return step (4.2), otherwise, export and treat reconstruct small echo high frequency coefficient matrix

(5) according to treating reconstruct small echo high frequency coefficient matrix calculate basic covariance matrix Ω and residual covariance matrix Π:

$\mathrm{\&Omega;}={\left(\mathrm{diag}0\left(\mathrm{Gamma}(a_0+c,a_0+\frac{1}{2}\mathrm{diag}({\tilde{X}}^T*\tilde{X}))\right)\right)}^{-1},$

$\mathrm{\&Pi;}={\left(\mathrm{diag}0\left(\mathrm{Gamma}(a_0+\frac{M}{2},a_0+\frac{1}{2}\mathrm{diag}\left({(Y-\mathrm{\&Phi;}*\tilde{X})}^T(Y-\mathrm{\&Phi;}*\tilde{X})\right))\right)\right)}^{-1},$

Wherein parameter a ₀for Q ties up constant vector, it is multiplication operations that the value of vectorial each element is 0.000001, *, ^tfor matrix transpose operation, ^-1for inversion operation, diag () operation represents takes out composition of vector by the diagonal entry of matrix in bracket, Gamma (,) operate expression generation one with Section 1 in bracket for form parameter, in bracket, Section 2 is the vector of the gamma distribution of scale parameter, and two in its bracket are Q with the vector produced and tie up, and diag0 () operation represents generation square formation, the diagonal entry of square formation is the vector element in bracket, and off diagonal element is 0;

(6) by outer iteration number of times t and model solution iteration total degree N ₁compare: if outer iteration number of times t≤N ₁, then t=t+1, returns step (4), otherwise performs step (7);

(7) previous generation's stack result is established equal total stack result O, by outer iteration number of times t respectively with model solution iteration total degree N ₁iteration total degree N is solved with coefficient ₂compare: if outer iteration times N ₁<t<N ₁+ N ₂, then total stack result is calculated outer iteration number of times t=t+1, returns step (4), otherwise calculates total stack result calculate final reconstruct small echo high frequency coefficient matrix perform step (8);

(8) according to the low frequency wavelet coefficient of dissociation L retained and final reconstruct small echo high frequency coefficient matrix carry out wavelet inverse transformation, obtain the restructuring graph of former figure.

The present invention catches the aggregation of small echo by setting up multivariate model Gauss, use based on the parameter adaptive correction of Gibbs sampling method to multivariate Gaussian model when solving, and instruct with the non-zero support of ambiguous location to wavelet coefficient that the wavelet low frequency coefficient retained is obtained by rim detection, the reconstruction quality of image and robustness are significantly improved.

Accompanying drawing explanation

Fig. 1 is general flow chart of the present invention;

Fig. 2 is the schematic diagram in the present invention, small echo high-frequency sub-band coefficient being converted into multivariate matrix of coefficients;

Fig. 3 is the schematic diagram determining fuzzy edge position in the present invention;

Fig. 4 be the present invention and prior art when sampling rate is 40% to the reconstruction result figure of Boat image;

Fig. 5 is that the Y-PSNR of the Lena image using the present invention and prior art to reconstruct is with sampling rate change curve.

Embodiment

Below in conjunction with the drawings and specific embodiments, the present invention is described in detail.

With reference to Fig. 1, specific embodiment of the invention step is as follows:

Step one, transmit leg sends observing matrix, low frequency wavelet coefficient of dissociation and high frequency multivariate calculation matrix.

(1a) image transmit leg is at wavelet field observed image, low frequency wavelet coefficient of dissociation is all retained as to the observation of low frequency wavelet coefficient of dissociation, with Random Orthogonal Gauss observing matrix Φ to horizontal high-frequent sub-band coefficients A ₁, vertical high frequency sub-band coefficients A ₂with diagonal angle high-frequency sub-band coefficient A ₃carry out multivariate observation respectively, obtain horizontal high-frequent subband multivariate calculation matrix Y ₁, vertical high frequency subband multivariate calculation matrix Y ₂with diagonal angle high-frequency sub-band multivariate calculation matrix Y ₃, the process of wherein multivariate observation is:

(1a1), by horizontal high-frequent sub-band coefficients A ₁, vertical high frequency sub-band coefficients A ₂, diagonal angle horizontal high-frequent sub-band coefficients A ₃m size is all divided into be individual block of pixels, each block of pixels pulls into row vector, and the M × Q changed into ties up horizontal multivariate matrix of coefficients X ₁, vertical multivariate matrix of coefficients X ₂, diagonal angle multivariate matrix of coefficients X ₃as shown in Figure 2, by the high-frequency sub-band coefficient A of 256 × 256 pixels according to 3 × 3 blocks of pixels, owing to can not the last two rows of whole point of high-frequency sub-band coefficient A and last two row replenish with zero, the block of pixels of every a line correspondence each 3 × 3 of multivariate matrix of coefficients X pulls into row;

(1a2), to horizontal multivariate matrix of coefficients X ₁, vertical multivariate matrix of coefficients X ₂, diagonal angle multivariate matrix of coefficients X ₃observe, obtain horizontal high-frequent subband multivariate calculation matrix Y ₁=Φ * X ₁, vertical high frequency subband multivariate calculation matrix Y ₂=Φ * X ₂, diagonal angle high-frequency sub-band multivariate calculation matrix Y ₃=Φ * X ₃;

(1b) image transmit leg sends Random Orthogonal Gauss observing matrix Φ, low frequency wavelet coefficient of dissociation L, horizontal high-frequent subband multivariate calculation matrix Y ₁, vertical high frequency subband multivariate calculation matrix Y ₂with diagonal angle high-frequency sub-band multivariate calculation matrix Y ₃to take over party.

Step 2, according to the observing matrix Φ, the low frequency wavelet coefficient of dissociation L that receive and high-frequency sub-band multivariate calculation matrix Y, obtains nonzero coefficient group index set u by the guidance of rim detection and correlativity.

(2a) take over party receives observing matrix Φ, low frequency wavelet coefficient of dissociation L, the horizontal high-frequent subband multivariate calculation matrix Y of image transmit leg transmission ₁, vertical high frequency subband multivariate calculation matrix Y ₂with diagonal angle high-frequency sub-band multivariate calculation matrix Y ₃, because the reconstruct mode of three high-frequency sub-band coefficients is the same, therefore three high-frequency sub-band multivariate calculation matrix unification Y represent;

(2b) high-frequency sub-band being all zero by the low frequency wavelet coefficient of dissociation L and three of reception carries out wavelet inverse transformation, obtains Edge obscuring image, as shown in the wavelet inverse transformation process in Fig. 3;

(2c) edge blurred picture Canny edge detection method carries out rim detection, obtains marginal position, and the pixel extracting position, corresponding edge in Edge obscuring image obtains fuzzy edge, as shown in the Canny edge detection process in Fig. 3;

(2d) one deck wavelet transformation is carried out to fuzzy edge, obtain fuzzy edge small echo high frequency coefficient and fuzzy edge wavelet low frequency coefficient, as shown in wavelet decomposition in Fig. 3; The position that fuzzy edge small echo high frequency coefficient absolute value is greater than threshold value h is set to 1, the position being less than threshold value h is set to 0 and obtains initial fuzzy location matrix, as obtained shown in ambiguous location matrix in Fig. 3, initial fuzzy location matrix is arranged in the Multivariable Fuzzy location matrix E of M × Q dimension according to the form of multi-variable matrix, wherein M is the columns of Φ, Q is the columns of Y, threshold value h=0.2;

(2e) to be multiplied with high-frequency sub-band multivariate calculation matrix Y the correlation matrix Φ obtained according to the transposition of observing matrix Φ ^t* Y, by the absolute value of this correlation matrix | Φ ^t* Y| and Multivariable Fuzzy location matrix E is weighted summation, obtain coefficient importance matrix V=| Φ ^t* Y|+w*E, wherein w is weighting coefficient;

(2f) every a line of coefficient importance matrix V is added, obtains coefficient importance vector, the index of c maximum in a coefficient importance vector element is taken out and forms nonzero coefficient group index set u={s ₁, s ₂..., s _i..., s _c, wherein i=1,2 ..., c, s _irepresent the index of i-th nonzero coefficient group, c is the integer being less than M.

Step 3, initialization iterations and multivariate Gaussian model.

(3a) initialization outer iteration number of times t=1, model solution iteration total degree is N ₁, it is N that coefficient solves iteration total degree ₂;

(3b) reconstruct small echo high frequency coefficient matrix is treated in initialization be the null matrix of M × Q dimension with total stack result O, wherein M is the columns of Φ, and Q is the columns of Y;

(3c) reconstruct small echo high frequency coefficient matrix is treated in initialization every a line all obey basic covariance matrix Ω multivariate Gaussian distribution, wherein Ω be Q × Q tie up matrix, its diagonal entry is 1, and non-is 0 to line element;

(3d) initialization residual error every a line all obey residual covariance matrix Π multivariate Gaussian distribution, wherein Π be Q × Q tie up matrix, its diagonal entry is 0.01, and off diagonal element is 0.

Step 4, according to observing matrix Φ, high-frequency sub-band multivariate calculation matrix Y, basic covariance matrix Ω, residual covariance matrix Π and nonzero coefficient group index set u, is calculated by Gibbs sampling method and treats reconstruct small echo high frequency coefficient matrix

(4a) initialization iterations i=1;

(4b) obtain treating reconstruct small echo high frequency coefficient matrix according to basic covariance matrix Ω, residual covariance matrix Π, observing matrix Φ s _ithe covariance matrix of row:

Wherein for the s of observing matrix Φ _irow, for transposition, ^-1for inversion operation;

(4c) according to residual covariance matrix Π, observing matrix Φ with treat reconstruct small echo high frequency coefficient matrix s _ithe covariance matrix of row obtain treating reconstruct small echo high frequency coefficient matrix s _ithe mean vector of row:

Wherein for the kth of observing matrix Φ arranges, x _kfor treating reconstruct small echo high frequency coefficient matrix row k, wherein k is less than or equal to M for being more than or equal to 1 and being not equal to s _iinteger;

(4d) according to treating reconstruct small echo high frequency coefficient matrix s _ithe mean vector of row and covariance matrix set up corresponding multivariate Gaussian model, calculate and treat reconstruct small echo high frequency coefficient matrix s _irow coefficient $x_{s_i}=\mathrm{Gaussian}({\mathrm{\&mu;}}_{s_i},{\mathrm{\&Lambda;}}_{s_i}),$ Wherein, represent that stochastic generation one is obeyed mean vector and is covariance matrix is multivariate Gaussian distribution vector;

(4e) iterations i and nonzero coefficient group index set u element number c is compared: if iterations i<c, then iterations i is from adding 1, return step (4.2), otherwise, export and treat reconstruct small echo high frequency coefficient matrix

Step 5, according to treating reconstruct small echo high frequency coefficient matrix calculate basic covariance matrix Ω and residual covariance matrix Π:

$\mathrm{\&Omega;}={\left(\mathrm{diag}0\left(\mathrm{Gamma}(a_0+c,a_0+\frac{1}{2}\mathrm{diag}({\tilde{X}}^T*\tilde{X}))\right)\right)}^{-1},$

$\mathrm{\&Pi;}={\left(\mathrm{diag}0\left(\mathrm{Gamma}(a_0+\frac{M}{2},a_0+\frac{1}{2}\mathrm{diag}\left({(Y-\mathrm{\&Phi;}*\tilde{X})}^T(Y-\mathrm{\&Phi;}*\tilde{X})\right))\right)\right)}^{-1},$

Wherein, parameter a ₀for Q ties up constant vector, it is multiplication operations that the value of vectorial each element is 0.000001, *, ^tfor matrix transpose operation, ^-1for inversion operation, diag () operation represents takes out composition of vector by the diagonal entry of matrix in bracket, Gamma (,) operate expression generation one with Section 1 in bracket for form parameter, in bracket, Section 2 is the vector of the gamma distribution of scale parameter, and two in its bracket are Q with the vector produced and tie up, and diag0 () operation represents generation square formation, the diagonal entry of square formation is the vector element in bracket, and off diagonal element is 0.

Step 6, by outer iteration number of times t and model solution iteration total degree N ₁compare: if outer iteration number of times t≤N ₁, then t=t+1, returns step 4, otherwise, perform step 7.

Step 7, if previous generation's stack result equal total stack result O, by outer iteration number of times t respectively with model solution iteration total degree N ₁iteration total degree N is solved with coefficient ₂compare: if outer iteration times N ₁<t<N ₁+ N ₂, then total stack result is calculated and establish outer iteration number of times t=t+1, return step 4; Otherwise, calculate total stack result calculate final reconstruct small echo high frequency coefficient matrix perform step 8.

Step 8, according to the low frequency wavelet coefficient of dissociation L retained and final reconstruct small echo high frequency coefficient matrix through wavelet inverse transformation, obtain reconstructed image.

Effect of the present invention can be further illustrated by following emulation.

1, simulated conditions: emulation of the present invention is at windows 7, SPI, CPU Intel (R) Core (TM) i5-3470, basic frequency 3.20GHz, software platform is that Matlab R2011b runs, emulation select be 512 × 512 four width standard testing natural image Lena, Peppers, Boat, Barbara, divide block size 32 × 32 for OMP algorithm, for MPA method and the inventive method use threshold value is the Canny rim detection of 0.15, c=1800, the absolute value of correlation matrix and weighted value ω=1.5 of ambiguous location matrix, iterations N ₁=100, N ₂=200.

2, content and result is emulated:

Emulation 1: fixed sample rate 40% time, by the present invention and existing OMP, MPA method, standard testing natural image is reconstructed under wavelet field, the reconstruct visual effect of image Boat as shown in Figure 4, wherein Fig. 4 (a) is the former figure of Boat, Fig. 4 (b) is the partial enlarged drawing of Fig. 4 (a), Fig. 4 (c), Fig. 4 (e) and Fig. 4 (g) are the restructuring graph of OMP, MPA and the inventive method respectively, and Fig. 4 (d), Fig. 4 (f) and Fig. 4 (h) are the partial enlarged drawing of Fig. 4 (c), Fig. 4 (e) and Fig. 4 (g) respectively.

As can be seen from restructuring graph and partial enlarged drawing, it is better that the marginal portion of reconstructed image of the present invention keeps, and the noise of smooth is also than the much less of the reconstructed image of OMP, MPA.

Fixed sample rate 40% time, carry out respectively reconstructing for five times to image Lena, Peppers, Barbara, Boat that four width sizes are 512*512 by the present invention and existing OMP, MPA method, the mean value of the Y-PSNR PSNR of five reconstruct reconstruction result is as shown in table 1.

Table 1 natural image size is the reconstruction result of 512*512, sampling rate 40%

As can be seen from Table 1, average PSNR value OMP, the MPA method of reconstructed image of the present invention are all high, show that the quality of reconstructed image is good.

Emulation 2: when sampling rate is respectively 30%, 35%, 40%, 45%, 50%, be being reconstructed of the Lena image of 512 × 512 to size under wavelet field by the present invention and existing OMP, BEPA method, the mean value of the Y-PSNR PSNR of 5 reconstruction result is as shown in table 2.

Table 2 Lena image uses the reconstruction result of OMP, MPA and the inventive method under different sampling rate

As can be seen from Table 2, the inventive method is all higher than at the Y-PSNR PSNR that sampling rate is the result figure that 30%, 35%, 40%, 45%, 50% time obtains the PSNR that OMP and MPA method obtains, and namely the reconstructed image quality of method of the present invention is higher than OMP and MPA method.

Obtain OMP according to table 2 data, as shown in Figure 5, the horizontal ordinate in Fig. 5 represents sampling rate to the trend map that the PSNR of the Lena image that MPA and the inventive method reconstruct changes with sampling rate, and ordinate represents Y-PSNR PSNR (dB) value.

As seen from Figure 5, the PSNR value of reconstruction result figure that obtains of the inventive method is apparently higher than additive method.

To sum up, the present invention can obtain reconstructed image clearly well, compared with other reconstructing methods existing, invention increases the reconstruction quality of image.

Claims (3)

1. under multivariate observation based on the image reconstructing method of edge constraint, comprise the steps:

(1) take over party receives Random Orthogonal Gauss observing matrix Φ, low frequency wavelet coefficient of dissociation L, the horizontal high-frequent subband multivariate calculation matrix Y of image transmit leg transmission ₁, vertical high frequency subband multivariate calculation matrix Y ₂with diagonal angle high-frequency sub-band multivariate calculation matrix Y ₃, three high-frequency sub-band multivariate calculation matrix unification Y are represented;

(2) according to the observing matrix Φ, the low frequency wavelet coefficient of dissociation L that receive and high-frequency sub-band multivariate calculation matrix Y, the set of nonzero coefficient group index is obtained by the guidance of rim detection and correlativity: u={s ₁, s ₂..., s _i..., s _c, wherein s _irepresent the index of i-th nonzero coefficient group, i=1,2 ..., c, c are the columns being less than Φ;

(3) initialization outer iteration number of times t=1, model solution iteration total degree is N ₁, it is N that coefficient solves iteration total degree ₂; Reconstruct small echo high frequency coefficient matrix is treated in initialization be the null matrix of M × Q dimension with total stack result O, wherein M is the columns of Φ, and Q is the columns of Y; Reconstruct small echo high frequency coefficient matrix is treated in initialization every a line all obey basic covariance matrix Ω multivariate Gaussian distribution; Initialization residual error every a line all obey residual covariance matrix Π multivariate Gaussian distribution;

(4) according to observing matrix Φ, high-frequency sub-band multivariate calculation matrix Y, basic covariance matrix Ω, residual covariance matrix Π and nonzero coefficient group index set u, calculated by Gibbs sampling method and treat reconstruct small echo high frequency coefficient matrix

(5) according to treating reconstruct small echo high frequency coefficient matrix calculate basic covariance matrix Ω and residual covariance matrix Π:

Wherein parameter a ₀for Q ties up constant vector, it is multiplication operations that the value of vectorial each element is 0.000001, *, ^tfor matrix transpose operation, ^-1for inversion operation, diag () operation represents takes out composition of vector by the diagonal entry of matrix in bracket, Gamma (,) operate expression generation one with Section 1 in bracket for form parameter, in bracket, Section 2 is the vector of the gamma distribution of scale parameter, and two in its bracket are Q with the vector produced and tie up, and diag0 () operation represents generation square formation, the diagonal entry of square formation is the vector element in bracket, and off diagonal element is 0;

(6) by outer iteration number of times t and model solution iteration total degree N ₁compare: if outer iteration number of times t≤N ₁, then t=t+1, returns step (4), otherwise performs step (7);

(7) previous generation's stack result is established equal total stack result O, by outer iteration number of times t respectively with model solution iteration total degree N ₁iteration total degree N is solved with coefficient ₂compare: if outer iteration times N ₁<t<N ₁+ N ₂, then total stack result is calculated outer iteration number of times t=t+1, returns step (4), otherwise calculates total stack result calculate final reconstruct small echo high frequency coefficient matrix perform step (8);

(8) according to the low frequency wavelet coefficient of dissociation L retained and final reconstruct small echo high frequency coefficient matrix carry out wavelet inverse transformation, obtain the restructuring graph of former figure.

2. according to claim 1 under multivariate observation based on the image reconstructing method of edge constraint, observing matrix Φ, low frequency wavelet coefficient of dissociation L according to reception described in step 2 and high-frequency sub-band multivariate calculation matrix Y, obtain nonzero coefficient group index set u by the guidance of rim detection and correlativity, concrete steps are as follows:

(2.1) high-frequency sub-band being all zero by the low frequency wavelet coefficient of dissociation L and three of reception carries out wavelet inverse transformation, obtains Edge obscuring image;

(2.2) edge blurred picture carries out rim detection, obtains marginal position;

(2.3) pixel extracting position, corresponding edge in Edge obscuring image obtains fuzzy edge;

(2.4) one deck wavelet transformation is carried out to fuzzy edge, obtain fuzzy edge small echo high frequency coefficient and fuzzy edge wavelet low frequency coefficient; The position that fuzzy edge small echo high frequency coefficient absolute value is greater than threshold value h is set to 1, the position being less than threshold value h is set to 0 and obtains initial fuzzy location matrix, initial fuzzy location matrix is arranged in the Multivariable Fuzzy location matrix E of M × Q dimension according to the form of multi-variable matrix, wherein M is the columns of Φ, Q is the columns of Y, threshold value h=0.2;

(2.5) to be multiplied with high-frequency sub-band multivariate calculation matrix Y the correlation matrix Φ obtained according to the transposition of observing matrix Φ ^t* Y, by the absolute value of this correlation matrix | Φ ^t* Y| and Multivariable Fuzzy location matrix E is weighted summation, obtain coefficient importance matrix V=| Φ ^t* Y|+w*E, wherein w is weighting coefficient;

(2.6) every a line of coefficient importance matrix V is added, obtains coefficient importance vector, the index of c maximum in a coefficient importance vector element is taken out and forms nonzero coefficient group index set u={s ₁, s ₂..., s _i..., s _c, wherein i=1,2 ..., c, s _irepresent the index of i-th nonzero coefficient group, c is the integer being less than M.

3. according to claim 1 under multivariate observation based on the image reconstructing method of edge constraint, described in step 4 according to observing matrix Φ, high-frequency sub-band multivariate calculation matrix Y, basic covariance matrix Ω, residual covariance matrix Π and nonzero coefficient group index set u, calculated by gibbs sampler and treat reconstruct small echo high frequency coefficient matrix concrete steps are as follows:

(4.1) initialization iterations i=1;

(4.2) obtain treating reconstruct small echo high frequency coefficient matrix according to basic covariance matrix Ω, residual covariance matrix Π, observing matrix Φ s _ithe covariance matrix of row:

Wherein for the s of observing matrix Φ _irow, for transposition

(4.3) according to residual covariance matrix Π, observing matrix Φ with treat reconstruct small echo high frequency coefficient matrix s _ithe covariance matrix Λ of row _si, obtain treating reconstruct small echo high frequency coefficient matrix s _ithe mean vector of row:

Wherein for the kth of observing matrix Φ arranges, x _kfor treating reconstruct small echo high frequency coefficient matrix row k, wherein k is less than or equal to M for being more than or equal to 1 and being not equal to s _iinteger;

(4.4) according to treating reconstruct small echo high frequency coefficient matrix s _ithe mean vector of row and covariance matrix set up corresponding multivariate Gaussian model, calculate and treat reconstruct small echo high frequency coefficient matrix s _irow coefficient wherein, represent that stochastic generation one is obeyed mean vector and is covariance matrix is multivariate Gaussian distribution vector;

(4.5) iterations i and nonzero coefficient group index set u element number c is compared: if iterations i<c, then iterations i is from adding 1, return step (4.2), otherwise, export and treat reconstruct small echo high frequency coefficient matrix