CN102183755A - Novel high-resolution orientation-estimating method based on Cauchy Gaussian model - Google Patents

Abstract

The invention relates to a novel high-resolution orientation-estimating method based on a Cauchy Gaussian model. By the method, the requirements of the underwater wideband high-frequency three-dimensional imaging of a buried mine detecting sonar on a target orientation resolving power are met. The method comprises the steps of: by using the Cauchy-Gaussian signal model, expressing an orientation spectrum as a regularized space-domain Fourier transform constraint optimizing problem; providing a spectral estimator tolerant to model parameters; by seeking sparse distribution features of a source signal in a space, realizing high-resolution spectral orientation estimation within a single-frequency-domain snapshot; and resolving coherent sources without decorrelation. The method has the advantages that: by the method, the high orientation resolution is achieved within the single snapshot through the constraint optimization of a cost function; the tolerance to parameters of a signal probability distributing model is high, so that the method is practical in engineering; the numerical simulation and sea trial data analyzing results prove that the method is applicable for the high target orientation resolution of a small-aperture subarray in shortage of the snapshot.

Description

A kind of novel high-resolution direction estimation method based on Cauchy-Gauss model

Technical field

The present invention relates to the sensor array signal process field, mainly is a kind of novel high-resolution direction estimation method based on Cauchy-Gauss model.

Background technology

The aspect estimation problem has caused extensive studies interest at sonar, radar and field such as communicate by letter all the time.Conventional wave beam forms and detects and the estimating target orientation by spatial matched filtering, but the target azimuth resolving power depends on the pore size of basic matrix wavelength unit.Adaptive beam forms [1] and subspace high-resolution direction of arrival (direction of arrival, DOA) algorithm for estimating as MUSIC, can obtain the high-resolution performance, but need enough signal to noise ratio (S/N ratio) conditions, and need long time average to improve the precision of covariance matrix.The signal subspace decomposition algorithm needs the pre-estimation signal number, when having coherent source, need separate relevant by the submatrix space smoothing.Thereby under restricted occasion in basic matrix aperture and high dynamic condition, conventional high resolution processing method is difficult to effective application.

The spatial domain constrained optimization is handled (spatial processing:optimized and constrained, SPOC) be a kind of novel single snap high resolution algorithm that occurred in recent years, the SPOC algorithm is based on the MAP criterion, constraint by signal model, optimal estimation goes out at each constantly, the sparse distribute power function of signal source in the space, thereby broken through the restriction of fast umber of beats and coherent source.People such as Korea S Kam are successfully applied to the seabed synthetic aperture sonar with the SPOC algorithm and lay the target direction high resolution processing.

The present invention utilizes wave beam formation and spatial domain Fourier transform at mathematical equivalence, studies a kind of iteration imparametrization azimuth spectrum estimation technique.Utilize Cauchy (Cauchy) to distribute and be similar to spatial domain DFT sample of signal, and signal minimum modulus constraint cost function is carried out regularization, the regularization cost function is optimized obtains spatial domain high-resolution DFT.Propose a kind of suboptimum, to the azimuth spectrum estimator of signal model parameters tolerance, to improve the robustness of algorithm.

Summary of the invention

The objective of the invention is to overcome the shortcoming and defect of prior art, a kind of novel high-resolution direction estimation method based on Cauchy-Gauss model is provided, bury of the requirement of thunder detection sonar the target azimuth resolving power to satisfy wideband high-frequency three-dimensional imaging under water.

The present invention solves the technical scheme that its technical matters adopts: this novel high-resolution direction estimation method based on Cauchy-Gauss model comprises the steps:

(1) basic matrix is received time-domain signal and transforms to frequency domain:

X(k)＝[X ₀(k)，...，X _M-1(k)] ^T(1)

N beamformer output signal of k frequency is expressed as spatial domain DFT form:

$\text{Y}(n,k)={\mathrm{\Σ}}_{m=0}^{M-1}{X}_m\left(k\right)e^{-j2\mathrm{\π}{f}_k\mathrm{md}\sin {\mathrm{\θ}}_n/c}$

$={\mathrm{\Σ}}_{m=0}^{N-1}{X}_m\left(k\right)e^{-j2\mathrm{\πmn}/N},N>M---\left(2\right)$

N counts for spatial domain DFT, and n discrete space orientation satisfies

θ _n＝asin{cn/(Nf _kd)}(3)

(2) utilize the character of DFT, space DFT spectrum estimated statement be shown as linear inversion find the solution form:

X(k)＝FY(k)(6)

Wherein, F is the Fourier transformation matrix, the beamformer output signal vector of Y (k) k frequency;

(3) structure canonical cost function:

$\tilde{J}={\mathrm{\Σ}}_{n=0}^{N-1}\log \left({\left|Y(n,k)\right|}^2\right)+{\left|\right|X-\mathrm{FY}\left|\right|}_2^2---\left(4\right)$

(4)Y＝(λQ ^-1+F ^HF) ^-1F ^HX(11)

Wherein, Q is a diagonal matrix:

$Q=\mathrm{diag}\{I_N+{\mathrm{YY}}^H/{\mathrm{\σ}}_Y^2\}---\left(12\right)$

Adjusting parameter lambda in the formula (11) is got fixing empirical value, again iterative high-resolution spatial spectrum in the following manner:

1) with Y ⁽⁰⁾=DFT{X} is as initial value, and substitution formula (17) obtains Q ⁽⁰⁾

$\tilde{Q}=\mathrm{diag}\left\{{\mathrm{YY}}^H\right\}---\left(17\right)$

2) the i time iteration has

Y ⁽ⁱ⁾＝Q ^(i-1)F ^H(λI _M+FQ ^(i-1)F ^H) ^-1X(5)

Utilize following formula to upgrade the Q matrix:

Q ⁽ⁱ⁾＝diag{Y ⁽ⁱ⁾Y ^(i)H}；

3) judge the condition of convergence:

|J ⁽ⁱ⁾-J ^(i-1)|/|J ^(i-1)|＜ε(6)

ε is an a small amount of, if do not satisfy condition, makes i=i+1, returns the 2nd) step; If satisfy condition, then export the spatial spectrum Y in current iteration cycle ⁽ⁱ⁾

(5) output is satisfied the high-resolution spatial spectrum of the condition of convergence to follow-up image display processing system.

The effect that the present invention is useful is: a kind of High Resolution DOA Estimation based on Cauchy-gaussian signal model that the present invention studied, realized orientation high-resolution in single snap by the constrained optimization of cost function.And this method is to the tolerance height of signal probability distributed model parameter, thereby has engineering practicability.Numerical simulation and sea examination data analysis presentation of results the inventive method are applicable to the target azimuth high-resolution of the small-bore basic matrix of owing under the snap condition.Similarly, the inventive method can be used for the high-resolution line spectrum estimation of naval vessels radiated noise.

Description of drawings

Fig. 1 is CBF F-K spectra (can't differentiate normalized frequency is two signals at 0.2 place);

Fig. 2 is RCG-SPEC algorithm F-K spectra (improved the robustness of algorithm, guaranteed the high-resolution performance simultaneously);

Embodiment

The invention will be further described below in conjunction with drawings and Examples:

This novel high-resolution direction estimation method based on Cauchy-Gauss model of the present invention may further comprise the steps:

Step 1: form is found the solution in the linear inversion that Estimation of Spatial Spectrum is expressed as spatial domain DFT:

X(k)＝FY(k)(6)

Wherein, X (k) is the basic matrix frequency-region signal vector of k frequency, and F is the Fourier transformation matrix, the beamformer output signal vector of Y (k) k frequency.

Step 2:, ask the constrained optimization regular solution of spatial domain DFT linear inversion: at first construct cost function based on Cauchy-gaussian signal model:

$J={\mathrm{\Σ}}_{n=0}^{N-1}\log (1+\frac{{\left|Y(n,k)\right|}^2}{{\mathrm{\σ}}_Y^2})+\frac{{\left|\right|X-\mathrm{FY}\left|\right|}_2^2}{{\mathrm{\σ}}_W^2}---\left(10\right)$

Wherein,

Be noise variance,

The variance of expression variable Y.To formula (10) differentiate, and to make it be zero, and the optimum solution that can get cost function is

Y＝(λQ ^-1+F ^HF) ^-1F ^HX(11)

Wherein,

Q is a diagonal matrix:

$Q=\mathrm{diag}\{I_N+{\mathrm{YY}}^H/{\mathrm{\σ}}_Y^2\}---\left(12\right)$

Step 3: model parameter toleranceization processing: regular terms in the formula (10) and minimum modulus item are got equal weight, and ignore constant term, write as a kind of canonical cost function form of suboptimum:

$\tilde{J}={\mathrm{\Σ}}_{n=0}^{N-1}\log \left({\left|Y(n,k)\right|}^2\right)+{\left|\right|X-\mathrm{FY}\left|\right|}_2^2---\left(16\right)$

As seen, the optimization solution procedure of formula (16) does not rely on the model parameter that Cauchy distributes, thereby has improved the robustness of algorithm.The Q matrix of formula (12) becomes:

$\tilde{Q}=\mathrm{diag}\left\{{\mathrm{YY}}^H\right\}---\left(17\right)$

Lambda parameter in the while hold mode (11) is to satisfy the orthotropicity of matrix.

Step 4: iterative high-resolution spatial spectrum:

The iterative step is as follows:

1) with Y ⁽⁰⁾=DFT{X} is as initial value, and substitution formula (17) obtains Q ⁽⁰⁾

2) the i time iteration has

Y ⁽ⁱ⁾＝Q ^(i-1)F ^H(λI _M+FQ ^(i-1)F ^H) ^-1X(18)

Utilize following formula to upgrade the Q matrix:

Q ⁽ⁱ⁾＝diag{Y ⁽ⁱ⁾Y ^(i)H}；

3) judge the condition of convergence:

|J ⁽ⁱ⁾-J ^(i-1)|/|J ^(i-1)|＜ε(19)

ε is an a small amount of.If do not satisfy condition, make i=i+1, return the 2nd) step; If satisfy condition, then export the spatial spectrum Y in current iteration cycle ⁽ⁱ⁾

Concrete steps are as follows:

Step 1: consider that array element distance is the M unit equal space line array of d, the frequency-region signal vector of k frequency can be expressed as:

X(k)＝[X ₀(k)，...，X _M-1(k)] ^T(7)

Wherein,

$X_m\left(k\right)={\mathrm{\Σ}}_{i=1}^{D}S\left(k\right)e^{j2\mathrm{\π}{f}_k\mathrm{md}\sin {\mathrm{\θ}}_i/c}+W_m\left(k\right)---\left(8\right)$

Be the signal of m array element at k frequency, m=0,1 ..., M-1, Wm (k) is the noise signal of corresponding k frequency, and c is the velocity of sound, and D is a signal number, and θ i is the orientation of i the relative beam of signal.

N beamformer output signal of k frequency can be expressed as

$\text{Y}(n,k)={\mathrm{\Σ}}_{m=0}^{M-1}{X}_m\left(k\right)e^{-j2\mathrm{\π}{f}_k\mathrm{md}\sin {\mathrm{\θ}}_n/c}$

$={\mathrm{\Σ}}_{m=0}^{N-1}{X}_m\left(k\right)e^{-j2\mathrm{\πmn}/N},N>M---\left(9\right)$

Following formula is a spatial domain zero padding DFT form, and N counts for spatial domain DFT, and n discrete space orientation satisfies

θ _n＝asin{cn/(Nf _kd)}(10)

Formula (3) and formula (4) show, by zero padding spatial domain DFT, can realize fast that the frequency domain wave beam forms.But similar to the conversion of frequency domain to time domain, by not improving bearing resolution in space dimension zero padding, the bearing resolution that conventional spectrum is estimated depends on effective physics pore size of array.

According to formula (3), its inverse transformation (IDFT) is

$X_m\left(k\right)=\frac{1}{N}{\mathrm{\Σ}}_{n=0}^{N-1}Y(n,k)e^{j2\mathrm{\πmn}/N}---\left(11\right)$

Be expressed as with linear equation system

X(k)＝FY(k)(12)

Wherein, F is the Fourier transformation matrix, its m, and n element is: [f] _{M, n}=(1/M) exp (j2 π mn/N).Formula (6) is promptly owed alignment system of equations problem, can find the solution by various linear inversion computing method.Promptly do not adopt conventional zero padding to handle yet, utilize limited physics array element data message, be finally inversed by the sparse distribute power of signal on dimensional orientation.Adopt suitable constrained optimization to handle, can improve bearing resolution, relative intensity that again can reflected signal.

Step 2:, ask the constrained optimization regular solution of spatial domain DFT linear inversion based on Cauchy-gaussian signal model:

Utilize the sparse distribution feature of space incident plane wave signal in the orientation is gathered, can carry out regularization, and can improve the resolution of dimensional orientation spectrum, complex magnitude spectrum in space is similar to regard as obeys Cauchy's distribution spatial domain DFT inverting:

$p\left\{Y(n,k)\right\}\∝\frac{1}{1+{\left|Y(n,k)\right|}^2/{\mathrm{\σ}}_Y^2}---\left(13\right)$

Correspondingly, the separate vectorial Y's of each element is distributed as p{Y (k) }=∏ _nP{Y (n, k) }.

The variance of expression variable Y.

Suppose that noise W obeys multiple Gaussian distribution The condition of basic matrix frequency-region signal is distributed as:

$p\left(X\right|Y,{\mathrm{\σ}}_W)=\frac{1}{{\left(\mathrm{\π}{\mathrm{\σ}}_W^2\right)}^M}{e}^{-{\left|\right|Y-\mathrm{FX}\left|\right|}_2^2/{\mathrm{\σ}}_W^2}---\left(14\right)$

According to Bayes (Bayes) criterion, the posteriority of vectorial Y distributes and can be expressed as:

$p\left(Y\right|X,{\mathrm{\σ}}_Y,{\mathrm{\σ}}_W)=\frac{p\left(Y\right|{\mathrm{\σ}}_Y)p\left(X\right|Y,{\mathrm{\σ}}_W)}{p(\text{X|Y,}{\mathrm{\σ}}_Y,{\mathrm{\σ}}_W)}---\left(15\right)$

Utilize formula (7) and formula (8), make the maximum maximum a posteriori probability of formula (9) separate (MAP) even also following formula cost function minimum:

$J={\mathrm{\Σ}}_{n=0}^{N-1}\log (1+\frac{{\left|Y(n,k)\right|}^2}{{\mathrm{\σ}}_Y^2})+\frac{{\left|\right|X-\mathrm{FY}\left|\right|}_2^2}{{\mathrm{\σ}}_W^2}---\left(16\right)$

Wherein, first is the regular terms according to Cauchy's distributed structure of angular spectrum amplitude, and as describing estimating of the sparse property of angular spectrum.According to

The degree of rarefication of size adjustment DFT inverting; Second is the minimum modulus error, also is the bound term of array signal itself, passes through noise power

Regulate the weight between regular terms and minimum modulus item.

To formula (10) differentiate, and to make it be zero, and the optimum solution that can get cost function is

Y＝(λQ ^-1+F ^HF) ^-1F ^HX(17)

Wherein,

Q is a diagonal matrix:

$Q=\mathrm{diag}\{I_N+{\mathrm{YY}}^H/{\mathrm{\σ}}_Y^2\}---\left(18\right)$

Formula (11) is that an amount of decrease minimum modulus is separated in form.Because Q is relevant with Y, thereby must be met the estimator of formula (11) by iterative Owing to have

F ^H(λI _M+FQF ^H)＝(λQ ^-1+F ^HF)QF ^H(19)

Utilize the orthotropicity of matrix in the bracket of both sides, have

(λQ ^-1+F ^HF) ^-1F ^H＝QF ^H(λI _M+FQF ^H) ^-1(20)

Thereby formula (11) can be written as

Y＝QF ^H(λI _M+FQF ^H) ^-1X(21)

Utilize the advantage of formula (15) replacement formula (11) to be only need find the solution M * M and tie up inverse of a matrix, reduced calculated amount.

Step 3: model parameter toleranceization processing:

Step 2 is based on the high-resolution Estimation of Spatial Spectrum algorithm of Cauchy-Gauss model, and its advantage is not need to utilize the multidata snap that covariance matrix is carried out smoothly, realizes the high-resolution DOA estimation in single-frequency numeric field data snap.But above-mentioned algorithm needs estimated parameter With

Value, its estimated accuracy will directly influence the azimuth spectrum estimation performance.When During much larger than angular spectrum big or small, then the output of the spectrum of formula (15) will deteriorate to conventional spatial spectrum, otherwise, suitable

Value can the sparse distribution feature of restoring signal in the space, thereby can tell the multiple goal orientation that is separated by very near more subtly.

According to above analysis, when

The time, often can improve resolving power, thereby utilize Regular terms in the formula (10) and minimum modulus item are got equal weight, and ignore constant term, the azimuth spectrum based on the tolerance of Cauchy-Gauss model that obtains a kind of suboptimum is estimated:

$\tilde{J}={\mathrm{\Σ}}_{n=0}^{N-1}\log \left({\left|Y(n,k)\right|}^2\right)+{\left|\right|X-\mathrm{FY}\left|\right|}_2^2---\left(22\right)$

As seen, the optimization solution procedure of formula (16) does not rely on the model parameter that Cauchy distributes, thereby has improved the robustness of algorithm.The Q matrix of formula (12) becomes:

$\tilde{Q}=\mathrm{diag}\left\{{\mathrm{YY}}^H\right\}---\left(23\right)$

And the λ in the hold mode (11), to satisfy the orthotropicity of matrix, utilizing probability distribution theory, λ can be similar to and get a fixing constant, the desirable 1/M of empirical value.

Step 4: iterative high-resolution spatial spectrum:

The iterative step is as follows:

1) with Y ⁽⁰⁾=DFT{X} (conventional wave beam formation) is as initial value, and substitution formula (17) obtains Q ⁽⁰⁾

2) the i time iteration has

Y ⁽ⁱ⁾＝Q ^(i-1)F ^H(λI _M+FQ ^(i-1)F ^H) ^-1X(24)

Utilize following formula to upgrade the Q matrix:

Q ⁽ⁱ⁾＝diag{Y ⁽ⁱ⁾Y ^(i)H}；

3) judge the condition of convergence:

|J ⁽ⁱ⁾-J ^(i-1)|/|J ^(i-1)|＜ε(25)

ε is an a small amount of.If do not satisfy condition, make i=i+1, return the 2nd) step; If satisfy condition, then export the spatial spectrum Y in current iteration cycle ⁽ⁱ⁾

Specific implementation of the present invention is:

(1) basic matrix is received time-domain signal and transforms to frequency domain:

X(k)＝[X ₀(k)，...，X _M-1(k)] ^T(26)

N beamformer output signal of k frequency is expressed as spatial domain DFT form:

$\text{Y}(n,k)={\mathrm{\Σ}}_{m=0}^{M-1}{X}_m\left(k\right)e^{-j2\mathrm{\π}{f}_k\mathrm{md}\sin {\mathrm{\θ}}_n/c}$

$={\mathrm{\Σ}}_{m=0}^{N-1}{X}_m\left(k\right)e^{-j2\mathrm{\πmn}/N},N>M---\left(27\right)$

N counts for spatial domain DFT, and n discrete space orientation satisfies

θ _n＝asin{cn/(Nf _kd)}(28)

(2) utilize the character of DFT, space DFT spectrum estimated statement be shown as linear inversion find the solution form:

X(k)＝FY(k)(6)

Wherein, F is the Fourier transformation matrix, the beamformer output signal vector of Y (k) k frequency.

(3) structure canonical cost function:

$\tilde{J}={\mathrm{\Σ}}_{n=0}^{N-1}\log \left({\left|Y(n,k)\right|}^2\right)+{\left|\right|X-\mathrm{FY}\left|\right|}_2^2---\left(29\right)$

(4) the adjusting parameter lambda in the formula (11) is got fixing empirical value (as 1/M), again iterative high-resolution spatial spectrum in the following manner:

1) with Y ⁽⁰⁾=DFT{X} (conventional wave beam formation) is as initial value, and substitution formula (17) obtains Q ⁽⁰⁾

2) the i time iteration has

Y ⁽ⁱ⁾＝Q ^(i-1)F ^H(λI _M+FQ ^(i-1)F ^H) ^-1X(30)

Utilize following formula to upgrade the Q matrix:

Q ⁽ⁱ⁾＝diag{Y ⁽ⁱ⁾Y ^(i)H}；

3) judge the condition of convergence:

|J ⁽ⁱ⁾-J ^(i-1)|/|J ^(i-1)|＜ε(31)

ε is an a small amount of.If do not satisfy condition, make i=i+1, return the 2nd) step; If satisfy condition, then export the spatial spectrum Y in current iteration cycle ⁽ⁱ⁾

(5) output is satisfied the high-resolution spatial spectrum of the condition of convergence to follow-up image display processing system.

The simulation analysis example:

Realistic model: 10 yuan of even linear arrays, signal are the sinusoidal plane waves of 3 unit amplitudes, space normalization wave number (fdsin θ/c) be respectively 0.2,0.25 and-0.25; Normalized frequency is 0.2,0.2 and 0.35.Signal time length is 100 samples, and the standard deviation of Gaussian noise is σ n=0.2, and every channel signal adds Hamming window and transforms to frequency domain.The F-K spectra that compares CBF and two kinds of algorithms of tolerance CG-SPEC algorithm (being designated as RCG-SPEC).

Fig. 1, Fig. 2 are the normalized frequency-wavenumber spectrum of two kinds of algorithms.As seen, it is two signals of 0.2 that CBF can't differentiate frequency, and secondary lobe is higher.As seen from Figure 2, RCG-SPEC algorithm proposed by the invention, although a small amount of clutter occurs, main lobe is widened to some extent, the robustness of algorithm is higher.

Sea examination data analysis:

Utilize the performance of the sea of high-resolution three-dimensional imaging under water examination data detection algorithm.Transmitting is linear frequency modulation, observes the 9 yuan of linear arrays of direction of vertically walking to navigate and composes at the arrowband at 12kHz frequency place echo bearing.CBF and RCG-SPEC algorithm are at the single snap degree of depth-orientation slice map of a certain position of walking to navigate, and bright spot is the bead that suspends in the water among the figure.As seen, institute of the present invention extracting method is compared with conventional wave beam formation method, has effectively improved bearing resolution.

Because the basic matrix aperture is subjected to the restriction of underwater towed-body size, and towed body is in the motion, conventional high resolution algorithm based on covariance matrix is difficult to effective application, and institute of the present invention extracting method need provide a kind of feasible solution route than the occasion of high resolution for physical pore size is limited.

In addition to the implementation, all employings are equal to the technical scheme of replacement or equivalent transformation formation, all drop on the protection domain of requirement of the present invention.

Claims (1)

1. the novel high-resolution direction estimation method based on Cauchy-Gauss model is characterized in that: comprise the steps:

(1) basic matrix is received time-domain signal and transforms to frequency domain:

X(k)＝[X ₀(k)，...，X _M-1(k)] ^T (1)

N beamformer output signal of k frequency is expressed as spatial domain DFT form:

$Y(n,k)={\mathrm{\Σ}}_{m=0}^{M-1}{X}_m\left(k\right)e^{-j2\mathrm{\π}{f}_k\mathrm{md}\sin {\mathrm{\θ}}_n/c}$

$={\mathrm{\Σ}}_{m=0}^{N-1}{X}_m\left(k\right)e^{-j2\mathrm{\πmn}/N},N>M---\left(2\right)$

N counts for spatial domain DFT, and n discrete space orientation satisfies

θ _n＝asin{cn/(Nf _kd)} (3)

(2) utilize the character of DFT, space DFT spectrum estimated statement be shown as linear inversion find the solution form:

X(k)＝FY(k) (6)

Wherein, F is the Fourier transformation matrix, the beamformer output signal vector of Y (k) k frequency;

(3) structure canonical cost function:

$\tilde{J}={\mathrm{\Σ}}_{n=0}^{N-1}\log \left({\left|Y(n,k)\right|}^2\right)+{\left|\right|X-\mathrm{FY}\left|\right|}_2^2---\left(4\right)$

(4)Y＝(λQ ^-1+F ^HF) ^-1F ^HX (11)

Wherein, $\mathrm{\λ}={\mathrm{\σ}}_W^2/{\mathrm{\σ}}_Y^2,$ Q is a diagonal matrix:

$Q=\mathrm{diag}\{I_N+{\mathrm{YY}}^H/{\mathrm{\σ}}_Y^2\}---\left(12\right)$

Adjusting parameter lambda in the formula (11) is got fixing empirical value, again iterative high-resolution spatial spectrum in the following manner:

1) with Y ⁽⁰⁾=DFT{X} is as initial value, and substitution formula (17) obtains Q ⁽⁰⁾

$\tilde{Q}=\mathrm{diag}\left\{{\mathrm{YY}}^H\right\}---\left(17\right)$

2) the i time iteration has

Y ⁽ⁱ⁾＝Q ^(i-1)F ^H(λI _M+FQ ^(i-1)F ^H) ^-1X (5)

Utilize following formula to upgrade the Q matrix:

Q ⁽ⁱ⁾＝diag{Y ⁽ⁱ⁾Y ^(i)H}；

3) judge the condition of convergence:

|J ⁽ⁱ⁾-J ^(i-1)|/|J ^(i-1)|＜ε (6)

ε is an a small amount of, if do not satisfy condition, makes i=i+1, returns the 2nd) step; If satisfy condition, then export the spatial spectrum Y in current iteration cycle ⁽ⁱ⁾

(5) output is satisfied the high-resolution spatial spectrum of the condition of convergence to follow-up image display processing system.

Images