METHOD FOR DIRECTION-OF-ARRIVAL ESTIMATION BASED ON SPARSE RECONSTRUCTION IN THE PRESENCE OF GAIN-PHASE ERROR

Abstract

Disclosed is a method for direction-of-arrival estimation based on sparse reconstruction in the presence of gain-phase error, which comprises the following steps: firstly, estimating a noise power and an gain error from an array received signal by adopting a characteristic decomposition method; then, based on a compensated covariance matrix, transforming a direction-of-arrival estimation problem into a non-convex optimization problem in a sparse frame by a method of sparse reconstruction; finally, estimating a grid angle and a deviation angle by using an alternate optimization method. This estimation method can effectively eliminate the influence of a phase error in direction-of-arrival estimation, and has better adaptability, which improves the resolution and estimation accuracy of the algorithm.

Claims

1. A method for direction-of-arrival estimation based on sparse reconstruction in the presence of gain-phase error, comprising the following steps: S1, calculating a covariance matrix R from an array received signal X(t), estimating a noise power by adopting a characteristic decomposition method, and estimating and compensating an gain error according to the noise power and main diagonal data of the covariance matrix to obtain a compensated covariance matrix R.sub.1; S2, according to the compensated covariance matrix R.sub.1 obtained in S1, transforming a direction-of-arrival estimation problem into a nonconvex optimization problem in a sparse frame by a method of sparse reconstruction, which is specifically realized through the following substeps: S2.1: according to the compensated covariance matrix R.sub.1, taking the magnitude of elements in the matrix to obtain |R.sub.1|, and taking the elements in an upper triangle area thereof, and eliminating the repeated elements of a same size in a main diagonal line, and then rearranging according to the following formula:
x=[|r.sub.1,1|,|r.sub.1,2|, . . . ,|r.sub.1,M|,|r.sub.2,3|, . . . ,|r.sub.2,M|, . . . ,|r.sub.M-1,M|].sup.T=|Bp|  (1)
B=[b(θ.sub.1),b(θ.sub.2), . . . ,b(θ.sub.K)]  (2)
p=[σ.sub.1.sup.2,σ.sub.2.sup.2, . . . ,σ.sub.K.sup.2].sup.T  (3) where B is a newly defined steering vector matrix composed of an angle θ.sub.k, p is a newly defined matrix composed of the power of K signals, σ.sub.k.sup.2 represents the power of a k.sup.th signal, (⋅).sup.T represents transposition, and b(θ.sub.k) represents a steering vector corresponding to the angle θ.sub.k, a value of which is shown in the following formula
b(θ.sub.k)=[1,e.sup.−j(τ.sup.k,2.sup.−τ.sup.k,1.sup.), . . . ,e.sup.−j(τ.sup.k,M.sup.−τ.sup.k,1.sup.),e.sup.−j(τ.sup.k,3.sup.−τ.sup.k,2.sup.), . . . ,e.sup.−j(τ.sup.k,M.sup.−τ.sup.k,2.sup.), . . . ,e.sup.−j(τ.sup.k,M.sup.−τ.sup.k,M-1.sup.)].sup.T  (4) where τ.sub.k,m represents a delay of the k.sup.th signal in an m.sup.th array element relative to a reference array element; S2.2: setting a space grid spacing Δ and constructing an overcomplete angle set Θ={−90°, −90°+Δ, . . . , 90° −Δ}, so as to extend the formula (1) to Θ to obtain an overcomplete output model of the following formula: $\begin{array}{cc}x=.Math.\mathrm{Bp}.Math.& \left(5\right)\\ B=\left[b\left(-90°\right),b\left(-90°+\Delta \right),\mathrm{.Math.},b\left(90°-\Delta \right)\right]& \left(6\right)\\ p=\{\begin{array}{cc}{\stackrel{\_}{p}}_k,& \theta ={\theta }_k\\ 0,& \mathrm{else}\end{array}& \left(7\right)\end{array}$ where B is a steering vector matrix formed by corresponding extension of B to Θ, and p is a matrix formed by corresponding extension of p to Θ; S2.3: if there is a deviation angle δ when an actual information source direction {tilde over (θ)} fails to fall strictly on the constructed grid, the first-order Taylor expansion is used to modify the steering vector B(θ) to
B({tilde over (θ)})=B(θ)+B′(θ).Math.δ  (8) where B({tilde over (ƒ)}) is the modified steering vector; S2.4: transforming the modified over-complete output model obtained in S2.3 into a nonconvex optimization problem of the following formula by an optimization theory
min.sub.p,δ∥x−|Bp+B′δp|∥.sub.2.sup.2∘  (9) S3, transforming a two-parameter non-convex optimization problem into a convex optimization problem by using an alternating optimization method, and obtaining a grid angle and a deviation angle by solving the convex optimization problem, and obtaining a final information source angle estimation value.

2. The method for direction-of-arrival estimation based on sparse reconstruction in the presence of gain-phase error according to claim 1, wherein S1 is implemented by the following sub step s: S1.1: calculating the covariance matrix R of the array received signal X(t), and then implementing eigenvalue decomposition on the covariance matrix R by using the following formula to obtain an eigenvalue λ.sub.m in a descending order:
R=Σ.sub.m=1.sup.Mλ.sub.mv.sub.mv.sub.m.sup.H  (10) where M represents a number of array elements, λ.sub.m represents the eigenvalue arranged in a descending order, v.sub.m represents an eigenvector corresponding to the eigenvalue λ.sub.m and (⋅).sup.H represents the conjugate transpose; S1.2: estimating the noise power {circumflex over (σ)}.sub.n.sup.2 by using the following formula according to the eigenvalue λ.sub.m obtained in S1.1, $\begin{array}{cc}{\hat{\sigma }}_n^2=\frac{1}{M-K}{.Math.}_{m=K+1}^M{\lambda }_m& \left(11\right)\end{array}$ where K represents a number of information sources; S1.3: estimating the gain error by using the following formula according to the obtained covariance matrix R and the estimated value {circumflex over (σ)}.sub.n.sup.2 of the noise power: $\begin{array}{cc}{\rho }_m=\sqrt{\frac{r_{m,m}-{\hat{\sigma }}_n^2}{r_{1,1}-{\hat{\sigma }}_n^2}}& \left(12\right)\end{array}$ where ρ.sub.m represents the estimated value of the gain error of the m.sup.th array element and r.sub.m,m represents the value at the covariance matrix (m, m); S1.4: compensating the estimated gain error matrix ρ.sub.m in the covariance matrix R by using the following formula, and eliminating the influence of the gain error to obtain a compensated covariance matrix R.sub.1:
R.sub.1=G.sup.−1(R−{circumflex over (σ)}.sub.n.sup.2I.sub.M)(G.sup.−1).sup.H  (13) where G=diag{[ρ.sub.1,ρ.sub.2, . . . , ρ.sub.M]} represents an gain error estimation matrix and I.sub.M represents an identity matrix with a size of M.

3. The method for direction-of-arrival estimation based on sparse reconstruction in the presence of gain-phase error according to claim 1, wherein S3 is implemented by the following substeps: S3.1: initializing a deviation angle matrix δ=0.sub.l, optimizing the problem of formula (13), and transforming the problem into the following formula:
min.sub.p,w∥w∥.sub.2.sup.2+γ.sub.1∥p∥.sub.2,1
s.t. p.sup.HA.sub.qp+w.sub.q=x.sub.q.sup.2  (14) wherein w=[w.sub.1, w.sub.2, . . . , w.sub.M].sup.T, γ.sub.1 represents a regularization constant, and, A.sub.q=b.sub.q.sup.Hb.sub.q, b.sub.q represents a q.sup.th line of B; S3.2: transforming formula (14) into a convex optimization problem of the following formula by using the idea of a feasible point pursuit algorithm, and solving formula (15) to obtain a sparse matrix p, and then obtaining the corresponding angle of a non-zero item in the sparse matrix p;
min.sub.p,w,c∥w∥.sub.2.sup.2+γ∥p∥.sub.2,1+μ.sub.1∥c∥.sub.1
s.t. p.sup.HA.sub.qp+w.sub.q≤x.sub.q.sup.2
2Re{z.sup.HA.sub.qp}+w.sub.q+c.sub.q≥x.sub.q.sup.2+z.sup.HA.sub.qz
p≥0
c.sub.q≥0  (15) where c=[c.sub.1, c.sub.2, . . . , c.sub.Q].sup.T, μ.sub.1 represents another regularization constant, and z represents an arbitrary matrix with the same specification with p; S3.3: solving the problem of formula (13) according to the sparse matrix p obtained in S3.2, and transforming the problem into the following problem: $\begin{array}{cc}{\min }_{\delta ,w}{.Math.w.Math.}_2^2+{\gamma }_2{.Math.\delta .Math.}_{2,1}s.t.C^2+2\mathrm{CD}\delta +{\delta }^H{E}_q\delta +w_q=x_q^2-\frac{\Delta }{2}\le \delta \le \frac{\Delta }{2}& \left(16\right)\end{array}$ where γ.sub.2 represents a regularization constant, C=Bp represents a known quantity, Dδ=B′δp, D represents an intermediate conversion quantity, δ represents a deviation angle matrix, and E.sub.q=d.sub.q.sup.Hd.sub.q, d.sub.q represents a q.sup.th line of D; S3.4: transforming the formula (16) into a convex optimization problem of the following formula by using the idea of a feasible point pursuit algorithm, and obtaining a deviation angle estimation matrix δ by solving the formula (17): $\begin{array}{cc}{\min }_{\delta ,w,c}{.Math.w.Math.}_2^2+{\gamma }_2{.Math.\delta .Math.}_{2,1}+{\mu }_2{.Math.c.Math.}_{1}s.t.C^2+2\mathrm{CD}\delta +{\delta }^H{E}_q\delta +w_q\le x_q^2{C}^2+2\mathrm{CD}\delta +c_q+2\mathrm{Re}\left\{z^H{E}_q\delta \right\}+w_q\ge x_q^2+z^H{E}_{q}z-\frac{\Delta }{2}\le 6\le \frac{\Delta }{2}{c}_q\ge 0& \left(17\right)\end{array}$ S3.5: obtaining an index matrix β corresponding to the grid angle matrix θ obtained in S3.2, and dot-multiplying a sum result of the grid angle matrix θ and the deviation angle matrix δ obtained in S3.4 with the index matrix β to obtain a final estimated source angle as follows:
{tilde over (θ)}=(θ+δ).Math.β  (18) where the index matrix β has a same dimension as the grid angle matrix θ, and the value of β at the index of the estimated angle is 1, with the rest being 0, (⋅) represents the dot multiplication of the matrix, that is, the multiplication of the corresponding elements of the matrix.

Description

BRIEF DESCRIPTION OF DRAWINGS

[0041] FIG. 1 is a flow chart of a method for DOA estimation based on sparse reconstruction in the presence of a gain-phase error.

[0042] FIG. 2 is a schematic diagram of grid division of an array spatial domain.

[0043] FIG. 3 is a comparison diagram of the relationship between the root mean square error and phase error in DOA estimation of the present disclosure and other algorithms in the same field.

[0044] FIG. 4 is a comparison chart of the relationship between the root mean square error and the signal-to-noise ratio in the DOA estimation of the present disclosure and other algorithms in the same field.

DESCRIPTION OF EMBODIMENTS

[0045] The purpose and effect of the present disclosure will become clearer from the following detailed description of the present disclosure according to the drawings and preferred embodiments. It should be understood that the specific embodiments described here are only used to explain, rather than to limit the present disclosure.

[0046] As shown in FIG. 1, the method for DOA estimation based on sparse reconstruction in the presence of a gain-phase error of the present disclosure includes the following steps:

[0047] S1, a covariance matrix is calculated from an array received signal, a noise power is estimated by adopting a characteristic decomposition method, and an gain error is estimated and compensated according to the noise power and main diagonal data of the covariance matrix to obtain a compensated covariance matrix; S1 is implemented by the following sub steps:

[0048] S1.1: calculating the covariance matrix R of the array received signal X(t), and then implementing eigenvalue decomposition on the covariance matrix R by using the following formula to obtain an eigenvalue λ.sub.m in a descending order

R=Σ.sub.m=1.sup.Mλ.sub.mv.sub.mv.sub.m.sup.H  (1)

[0049] where M represents a number of array elements, λ.sub.m represents the eigenvalue arranged in a descending order, v.sub.m represents an eigenvector corresponding to the eigenvalue λ.sub.m and (⋅).sup.H represents the conjugate transpose;

[0050] S1.2: estimating the noise power {circumflex over (σ)}.sub.n.sup.2 by using the following formula according to the eigenvalue λ.sub.m obtained in S1.1,

[00006]$\begin{array}{cc}{\hat{\sigma }}_n^2=\frac{1}{M-K}{.Math.}_{m=K+1}^M{\lambda }_m& \left(2\right)\end{array}$

[0051] where K represents a number of information sources;

[0052] S1.3: estimating the gain error by using the following formula according to the obtained covariance matrix R and the estimated value {circumflex over (σ)}.sub.n.sup.2 of the noise power

[00007]$\begin{array}{cc}{\rho }_m=\sqrt{\frac{r_{m,m}-{\hat{\sigma }}_n^2}{r_{1,1}-{\hat{\sigma }}_n^2}}& \left(3\right)\end{array}$

[0053] where ρ.sub.m represents the estimated value of the gain error of the mth array element and r.sub.m,m represents the value at the covariance matrix (m, m);

[0054] S1.4: compensating the estimated gain error matrix ρ.sub.m in the covariance matrix R by using the following formula, and eliminating the influence of the gain error to obtain a compensated covariance matrix R.sub.1

R.sub.1=G.sup.−1(R−{circumflex over (σ)}.sub.n.sup.2I.sub.M)(G.sup.−1).sup.H  (4)

[0055] where G=diag{[ρ.sub.1, ρ.sub.2, . . . , ρ.sub.M]} represents an gain error estimation matrix and I.sub.M represents an identity matrix with a size of M.

[0056] S2, according to the compensated covariance matrix obtained in S1, a direction-of-arrival estimation problem is transformed into a nonconvex optimization problem in a sparse frame by a method of sparse reconstruction; S2 is specifically realized through the following substeps:

[0057] S2.1: according to the compensated covariance matrix R.sub.1, taking the magnitude of elements in the matrix to obtain |R.sub.1|, and taking the elements in an upper triangle area thereof, and eliminating the repeated elements of a same size in a main diagonal line, and then rearranging according to the following formula

x=[|r.sub.1,1|,|r.sub.1,2|, . . . ,|r.sub.1,M|,|r.sub.2,3|, . . . ,|r.sub.2,M|, . . . ,|r.sub.M-1,M|].sup.T=|Bp|  (5)

B=[b(θ.sub.1),b(θ.sub.2), . . . ,b(θ.sub.K)]  (6)

p=[σ.sub.1.sup.2,σ.sub.2.sup.2, . . . ,σ.sub.K.sup.2].sup.T  (7)

[0058] where B is a newly defined steering vector matrix composed of an angle θ.sub.k, p is a newly defined matrix composed of the power of K signals, σ.sub.k.sup.2 represents the power of a k th signal, (⋅).sup.T represents transposition, and b(θ.sub.k) represents a steering vector corresponding to the angle θ.sub.k, a value of which is shown in the following formula

b(θ.sub.k)=[1,e.sup.−j(τ.sup.k,2.sup.−τ.sup.k,1.sup.), . . . ,e.sup.−j(τ.sup.k,M.sup.−τ.sup.k,1.sup.),e.sup.−j(τ.sup.k,3.sup.−τ.sup.k,2.sup.), . . . ,e.sup.−j(τ.sup.k,M.sup.−τ.sup.k,2.sup.), . . . ,e.sup.−j(τ.sup.k,M.sup.−τ.sup.k,M-1.sup.)].sup.T  (8)

[0059] where τ.sub.k,m represents a delay of the kth signal in a mth array element relative to a reference array element;

[0060] S2.2: setting a space grid spacing Δ and constructing an overcomplete angle set Θ={−90°, −90°+Δ, . . . , 90° −Δ}, so as to extend the formula (1) to Θ to obtain an overcomplete output model of the following formula

[00008]$\begin{array}{cc}x=.Math.\mathrm{Bp}.Math.& \left(9\right)\\ B=\left[b\left(-90°\right),b\left(-90°+\Delta \right),\mathrm{.Math.},b\left(90°-\Delta \right)\right]& \left(10\right)\\ p=\{\begin{array}{cc}{\stackrel{\_}{p}}_k,& \theta ={\theta }_k\\ 0,& \mathrm{else}\end{array}& \left(11\right)\end{array}$

[0061] where B is a steering vector matrix formed by corresponding extension of B to Θ, and p is a matrix formed by corresponding extension of p, to Θ;

[0062] S2.3: if there is a deviation angle δ when an actual information source direction {tilde over (θ)} fails to fall strictly on the constructed grid, the first-order Taylor expansion is used to modify the steering vector B(θ) to

B({tilde over (θ)})=B(θ)+B′(θ).Math.δ  (12)

[0063] where B({tilde over (θ)}) is the modified steering vector;

[0064] S2.4: transforming the modified over-complete output model obtained in S2.3 into a nonconvex optimization problem of the following formula by an optimization theory

min.sub.p,δ∥x−|Bp+B′δp|∥.sub.2.sup.2∘  (13).

[0065] S3, a two-parameter non-convex optimization problem is transformed into a convex optimization problem by using an alternating optimization method, and obtaining a grid angle and a deviation angle by solving the convex optimization problem, and obtaining a final information source angle estimation value; S3 is implemented by the following substeps:

[0066] S3.1: initializing a deviation angle matrix δ=0.sub.l, optimizing the problem of formula (13), and transforming the problem into the following formula

min.sub.p,w∥w∥.sub.2.sup.2+γ.sub.1∥p∥.sub.2,1

s.t. p.sup.HA.sub.qp+w.sub.q=x.sub.q.sup.2  (14)

[0067] wherein w=[w.sub.1, w.sub.2, . . . w.sub.M].sup.T, γ.sub.1 represents a regularization constant, and, A.sub.q=b.sub.q.sup.Hb.sub.q, b.sub.q represents a qth line of B;

[0068] S3.2: transforming formula (14) into a convex optimization problem of the following formula by using the idea of a feasible point pursuit algorithm, and solving formula (15) to obtain a sparse matrix p, and then obtaining the corresponding angle of a non-zero item in the sparse matrix p;

min.sub.p,w,c∥w∥.sub.2.sup.2+γ∥p∥.sub.2,1+μ.sub.1∥c∥.sub.1

s.t. p.sup.HA.sub.qp+w.sub.q≤x.sub.q.sup.2

2Re{z.sup.HA.sub.qp}+w.sub.q+c.sub.q≥x.sub.q.sup.2+z.sup.HA.sub.qz

p≥0

c.sub.q≥0  (15)

[0069] where c=[c.sub.1, c.sub.2, . . . , c.sub.Q].sup.T, μ.sub.1 represents another regularization constant, and z represents an arbitrary matrix with the same specification with p;

[0070] S3.3: solving the problem of formula (13) according to the sparse matrix p obtained in S3.2, and transforming the problem into the following problem

[00009]$\begin{array}{cc}{\min }_{\delta ,w}{.Math.w.Math.}_2^2+{\gamma }_2{.Math.\delta .Math.}_{2,1}s.t.C^2+2\mathrm{CD}\delta +{\delta }^H{E}_q\delta +w_q=x_q^2-\frac{\Delta }{2}\le \delta \le \frac{\Delta }{2}& \left(16\right)\end{array}$ [0071] where γ.sub.2 represents a regularization constant, C=Bp represents a known quantity, Dδ=B′δp, D represents an intermediate conversion quantity, δ represents a deviation angle matrix, and E.sub.q=d.sub.g.sup.Hd.sub.q, d.sub.q represents a qth line of D;

[0072] S3.4: transforming the formula (16) into a convex optimization problem of the following formula by using the idea of a feasible point pursuit algorithm, and obtaining a deviation angle estimation matrix δ by solving the formula (17)

[00010]$\begin{array}{cc}{\min }_{\delta ,w,c}{.Math.w.Math.}_2^2+{\gamma }_2{.Math.\delta .Math.}_{2,1}+{\mu }_2{.Math.c.Math.}_{1}s.t.C^2+2\mathrm{CD}\delta +{\delta }^H{E}_q\delta +w_q\le x_q^2{C}^2+2\mathrm{CD}\delta +c_q+2\mathrm{Re}\left\{z^H{E}_q\delta \right\}+w_q\ge x_q^2+z^H{E}_{q}z-\frac{\Delta }{2}\le 6\le \frac{\Delta }{2}{c}_q\ge 0& \left(17\right)\end{array}$

[0073] S3.5: obtaining an index matrix β corresponding to the grid angle matrix θ obtained in S3.2, and dot-multiplying a sum result of the grid angle matrix θ and the deviation angle matrix δ obtained in S3.4 with the index matrix β to obtain a final estimated source angle as follows

{tilde over (θ)}=(θ+δ).Math.β  (18)

[0074] where the index matrix β has a same dimension as the grid angle matrix θ, and the value of β at the index of the estimated angle is 1, with the rest being 0, (⋅) represents the dot multiplication of the matrix, that is, the multiplication of the corresponding elements of the matrix.

[0075] FIG. 2 is a schematic diagram of grid division of an array spatial domain, in which diamonds represent array elements, open circles represent grid points dividing the spatial domain, with a grid spacing being Δ, and filled circles represent actual directions of signals. When the hollow circle coincides with the solid circle, it means that the actual direction of the signal just falls on the grid, otherwise, the grid division model will produce a certain deviation error δ.

[0076] FIG. 3 is a comparison diagram of the relationship between the root mean square error and phase error in DOA estimation of the present disclosure and other algorithms in the same field. It can be seen from FIG. 3 that with the increase of an initial phase error, the root mean square error in DOA estimation of the present disclosure does not change, and this method (the proposed curve in the figure) can effectively eliminate the influence of a phase error in DOA estimation.

[0077] FIG. 4 is a comparison chart of the relationship between the root mean square error and the signal-to-noise ratio in DOA estimation between the present disclosure and other algorithms in the same field. It can be seen from FIG. 4 that the root mean square error of DOA estimation decreases with the increase of the signal-to-noise ratio, especially when the signal-to-noise ratio is greater than 15 dB, and the root mean square error of this method (the proposed curve in the figure) is smaller as compared with other algorithms, which shows that this method can improve the accuracy of DOA estimation.

[0078] It can be understood by those skilled in the art that the above description is only the preferred examples of the present disclosure, and is not intended to limit the present disclosure. Although the present disclosure has been described in detail with reference to the foregoing examples, those skilled in the art can still modify the technical solutions described in the foregoing examples or replace some of their technical features equivalently. Within the spirit and principle of the present disclosure, the modifications, equivalent replacements and so on shall be included within the scope of protection of the present disclosure.