CN113534040A - Coherent source-isolated gate DOA estimation method based on weighted second-order sparse Bayes - Google Patents

Coherent source-isolated gate DOA estimation method based on weighted second-order sparse Bayes Download PDF

Info

Publication number
CN113534040A
CN113534040A CN202110599238.0A CN202110599238A CN113534040A CN 113534040 A CN113534040 A CN 113534040A CN 202110599238 A CN202110599238 A CN 202110599238A CN 113534040 A CN113534040 A CN 113534040A
Authority
CN
China
Prior art keywords
signal
matrix
unitary transformation
algorithm
alpha
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN202110599238.0A
Other languages
Chinese (zh)
Other versions
CN113534040B (en
Inventor
沈明威
胥翔竣
蒋意扬
邱存银
俞帆
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Hohai University HHU
Original Assignee
Hohai University HHU
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Hohai University HHU filed Critical Hohai University HHU
Priority to CN202110599238.0A priority Critical patent/CN113534040B/en
Publication of CN113534040A publication Critical patent/CN113534040A/en
Application granted granted Critical
Publication of CN113534040B publication Critical patent/CN113534040B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S3/00—Direction-finders for determining the direction from which infrasonic, sonic, ultrasonic, or electromagnetic waves, or particle emission, not having a directional significance, are being received
    • G01S3/02—Direction-finders for determining the direction from which infrasonic, sonic, ultrasonic, or electromagnetic waves, or particle emission, not having a directional significance, are being received using radio waves
    • G01S3/14—Systems for determining direction or deviation from predetermined direction
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S7/00—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00
    • G01S7/02—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S13/00
    • G01S7/41—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S13/00 using analysis of echo signal for target characterisation; Target signature; Target cross-section
    • G01S7/418—Theoretical aspects
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F18/00—Pattern recognition
    • G06F18/20—Analysing
    • G06F18/23—Clustering techniques
    • G06F18/232—Non-hierarchical techniques
    • G06F18/2321—Non-hierarchical techniques using statistics or function optimisation, e.g. modelling of probability density functions
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F18/00—Pattern recognition
    • G06F18/20—Analysing
    • G06F18/24—Classification techniques
    • G06F18/241—Classification techniques relating to the classification model, e.g. parametric or non-parametric approaches
    • G06F18/2415—Classification techniques relating to the classification model, e.g. parametric or non-parametric approaches based on parametric or probabilistic models, e.g. based on likelihood ratio or false acceptance rate versus a false rejection rate
    • G06F18/24155—Bayesian classification
    • Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02D—CLIMATE CHANGE MITIGATION TECHNOLOGIES IN INFORMATION AND COMMUNICATION TECHNOLOGIES [ICT], I.E. INFORMATION AND COMMUNICATION TECHNOLOGIES AIMING AT THE REDUCTION OF THEIR OWN ENERGY USE
    • Y02D30/00—Reducing energy consumption in communication networks
    • Y02D30/70—Reducing energy consumption in communication networks in wireless communication networks

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Data Mining & Analysis (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Evolutionary Biology (AREA)
  • Evolutionary Computation (AREA)
  • General Engineering & Computer Science (AREA)
  • Artificial Intelligence (AREA)
  • Probability & Statistics with Applications (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Remote Sensing (AREA)
  • Computer Networks & Wireless Communication (AREA)
  • Measurement Of Velocity Or Position Using Acoustic Or Ultrasonic Waves (AREA)

Abstract

The invention discloses a coherent source outlier DOA estimation method based on weighted second-order sparse Bayes, which comprises the following steps: step S1: for received signal YcUnitary transformation realizes real number domain transformation and signal decorrelation processing, thereby obtaining a receiving signal Y after unitary transformation; step S2: performing singular value decomposition on the received signal Y after unitary transformation to obtain YsvAnd SsvCombining the MUSIC algorithm to construct a weighting vector w; step S3: solving hyper-parameter alpha for parameter matrix unitary transformation of overcomplete dictionary constructed by second-order guide vector and combining with weighted sparse Bayesian inference (OGWSBI)0Alpha, and offset beta, and outputs a DOA estimate
Figure DDA0003092317400000011
Simulation (Emulation)The experimental results show that: compared with the OGWSBI algorithm, under a coherent information source, the method provided by the invention can improve the estimation precision by about 50% and improve the operation efficiency by about 60%.

Description

Coherent source-isolated gate DOA estimation method based on weighted second-order sparse Bayes
Technical Field
The invention belongs to the field of radar signal processing, and particularly relates to a coherent source outlier DOA estimation method based on weighted second-order sparse Bayes.
Background
Direction Of Arrival (DOA) estimation is a fundamental problem Of array signal processing, and is one Of the important tasks in many fields such as radar and sonar. In the development process of the DOA estimation technology, a plurality of algorithms are proposed in sequence, including a traditional DOA estimation algorithm (Capon algorithm, forward prediction algorithm, maximum entropy algorithm, minimum modulus algorithm, and the like), a subspace-based algorithm, maximum likelihood estimation, and the like. The estimation accuracy of the subspace algorithm completely depends on whether the signal subspace and the noise subspace can be correctly represented, and the estimation performance of the method is greatly reduced under the conditions of reduced snapshot number, reduced signal-to-noise ratio, correlation between target signals and the like.
With the rise of the Compressed Sensing (CS) theory, researchers apply the sparse representation theory to DOA estimation and propose many algorithms for the sparsity of the spatial domain. Among them, the most representative algorithm is the L1-SVD algorithm. However, the above algorithm is based on the assumption that the target angle of the signal source is exactly on the observation matrix, and the estimation performance is degraded when there is a deviation, i.e., an off-grid problem.
In order to reduce the error caused by grid offset, researchers have conducted a lot of research on the off-grid scene, and in recent years, a DOA estimation method for the off-grid problem has been proposed. The off-grid sparse Bayesian inference (OGSBI) algorithm proposed by Yang et al in the literature based on Bayesian thought is representative, and the algorithm realizes high-precision DOA estimation under the coarse grid condition. Zhang et al, on the basis of OGSBI algorithm, combines MUSIC algorithm to construct prior weight vector of information source, and proposes weighted sparse Bayes (OGWSBI) algorithm. Compared with the OGSBI algorithm, the OGWSBI algorithm improves the convergence rate and the estimation precision during iteration. However, when the OGWSBI algorithm estimates the coherent signal source, due to the limitation of the MUSIC algorithm, the DOA estimation is missed in the spatial spectrum, so that the position where the signal source exists cannot be effectively weighted, and the estimation accuracy of the parameters is reduced. Therefore, new weighting coefficients need to be constructed to provide a correct a priori information for the iteration. In addition, due to the fact that high-order terms are absent in the Taylor expansion of the guide vector, the parameter estimation precision is reduced, and therefore the high-order terms need to be added, and the estimation precision is further improved.
Disclosure of Invention
The invention provides a coherent information source grid DOA estimation method based on weighted second-order sparse Bayes, aiming at solving the problems that high-order terms in guiding vector Taylor expansion are lacked and the DOA estimation precision is not high when the weighted sparse Bayes algorithm is applied to a coherent information source. On the basis of an OGWSBI algorithm, firstly, in order to reduce a model error caused by first-order Taylor expansion of a guide vector, the guide vector is expanded to second-order Taylor expansion; and then, unitary transformation is utilized to perform de-coherence processing on the received signals, recover the rank of the covariance matrix, and convert the estimation model from a complex number domain to a real number domain, thereby improving the estimation precision and the operation efficiency of the algorithm.
In order to achieve the purpose, the invention adopts the following technical scheme:
a coherent source outlier (DOA) estimation method based on weighted second-order sparse Bayes is used for the following models: receive a signal of
Figure BDA0003092317380000021
The sparse model expression is YcPhi (beta) S + N, wherein
Figure BDA0003092317380000022
Is a signal source, beta is an angular deviation,
Figure BDA0003092317380000023
representing noise, wherein M is an array element number, N is a grid division number, T is a fast beat number, K is a signal source number, and phi is an over-complete dictionary formed by guide vectors;
the DOA estimation method comprises the following steps:
step S1: for received signal YcUnitary transformation real number domain conversion and signal de-coherent processing are carried out, so that a receiving signal Y after unitary transformation is obtained;
step S2: decomposing the singular value of Y of the received signal after unitary transformation to obtain YsvAnd SsvCombining the MUSIC algorithm to construct a weighting vector w;
step S3: solving the hyperparameter alpha by combining the weighted sparse Bayesian inference to the unitary transformation of the parameter matrix of the overcomplete dictionary constructed by the second-order guide vector0Alpha, and offset beta, and outputs a DOA estimate
Figure BDA0003092317380000024
Further, step S3 includes the steps of:
step S31: constructing an over-complete dictionary by using a real source arrival angle guide vector of second-order Taylor expansion:
Figure BDA0003092317380000025
wherein ,AcFor the observation matrix:
Figure BDA0003092317380000026
Figure BDA0003092317380000027
is distributed in [0, π]Discrete sampling values of the space angles are obtained, and N is the number of grid divisions; the guide vector of the second-order Taylor expansion is
Figure BDA0003092317380000028
Figure BDA0003092317380000029
At an angle theta to the true signalkThe closest value of the grid to the grid value,
Figure BDA00030923173800000210
a (theta) is in
Figure BDA00030923173800000211
To the first order ofThe derivative(s) of the signal(s),
Figure BDA00030923173800000212
a (theta) is in
Figure BDA00030923173800000213
Second derivative of (a) of
Figure BDA00030923173800000214
The matrix is composed of
Figure BDA00030923173800000215
Figure BDA00030923173800000216
Composed matrix
Figure BDA00030923173800000217
Step S32: for A in overcomplete dictionaryc、BcAnd CcOver-complete dictionary after unitary transformation is obtained by performing unitary transformation
Figure BDA0003092317380000031
Step S33: by using
Figure BDA0003092317380000032
Λ＝diag(α)、p(αn)＝Γ(α|1,wnRho) which is a user selectivity parameter and a signal source SsvA posterior probability distribution of
Figure BDA0003092317380000033
Figure BDA0003092317380000034
Obtaining the mean value and the variance mu (t) ═ alpha of the signal0∑ΦHysv(t),t＝1,....K，∑＝(α0ΦHΦ+Λ-1)-1Updating the mean value mu (t) and the variance sigma of the signal with the overcomplete dictionary phi after the unitary transformation and based on the mean value mu (t) and the variance sigma of the signalVariance sigma updating hyperparameter alpha0α and offset β; judging whether the maximum iteration number i is reachedmax or ||αi+1-αi||2/||αi 2If the two are not satisfied, updating the overcomplete dictionary phi by using the updated offset beta, and continuously updating the mean value mu (t) and the variance sigma of the signal by using the updated overcomplete dictionary phi so as to update the overcompensate parameter alpha0α and offset β; if one or two conditions are met, the iteration is ended, and the DOA estimated value is output
Figure BDA0003092317380000035
Further, the updated hyper-parameter α0And alpha is respectively:
Figure BDA0003092317380000036
Figure BDA0003092317380000037
wherein, μ ═ μ (1),.. μ (K)]＝α0∑ΦHYsv，
Figure BDA0003092317380000038
ρ＝ρ/K，
Figure BDA0003092317380000039
b,c→0。
Further, the updated offset β is:
Figure BDA00030923173800000310
wherein ,
Figure BDA00030923173800000311
Figure BDA00030923173800000312
Figure BDA0003092317380000041
further, step S1 includes the following steps:
step S11: for received signal YcConstructing an augmentation matrix to form a central Hermite matrix:
Figure BDA0003092317380000042
wherein ,Yc *Is YcThe complex conjugate of (a) and (b),
Figure BDA0003092317380000043
step S12: centering Hermite matrix YaugUnitary transformation processing:
Figure BDA0003092317380000044
wherein Y represents YaugAfter unitary transformation, transforming the complex number domain into a matrix of a real number domain, wherein U is a unitary matrix;
when M is even, UMCan be expressed as:
Figure BDA0003092317380000045
wherein, I is an M/2 xM/2 unit matrix, J is an M/2 xM/2 exchange matrix, the minor diagonal of the exchange matrix is 1, and other elements are 0;
when M is odd, UMCan be expressed as:
Figure BDA0003092317380000046
wherein the dimensions of I and J are (M-1)/2;
in the same way, U2TIs a 2T-dimensional unitary matrix.
Further, in step S32, A in the overcomplete dictionary is checkedc、BcAnd CcA, B and C for unitary transformation are shown below:
Figure BDA0003092317380000047
Figure BDA0003092317380000048
Figure BDA0003092317380000049
further, step S2 includes the steps of:
performing singular value decomposition on the received signal Y after unitary transformation, namely Y ═ Us∑sVs T+Ue∑eVe TTo obtain a signal subspace UsNoise subspace Ue(ii) a Order to
Figure BDA00030923173800000410
Figure BDA00030923173800000411
For the real signal source, Y is obtainedsv、Ssv；
The inverse MUSIC spatial spectrum formula is:
Figure BDA0003092317380000051
calculate the weight value w ═ w1,…,wN]TMeanwhile, normalizing the weight to obtain a weighted vector w: w ═ w/min (w).
Has the advantages that: performing coherent demodulation processing on the received signals by using unitary transformation, and recovering the rank of the covariance matrix; meanwhile, in order to reduce the model error caused by the first-order Taylor expansion of the guide vector, the guide vector is expanded to the second-order Taylor expansion, and an over-complete word formed by the guide vector is formedUnitary transformation is carried out to realize a sparse model Y of the received signalcAnd real quantization of phi (beta) S + N, and finally solving the offset by combining with weighted sparse Bayes inference.
Drawings
FIG. 1 is a flow chart of a weighted second-order sparse Bayesian-based coherent source-leaving-fence DOA estimation method (IOGWSBI) in an embodiment of the present invention;
FIG. 2 is a cross-sectional view of angle estimation of the OGWSBI algorithm under the condition of a first-order Taylor expansion model of a steering vector under coherent and incoherent sources;
FIG. 3 is a cross-sectional view of coherent source angle estimation of MUSIC algorithm before and after unitary transformation processing;
FIG. 4 is a cross-sectional view of coherent source angle estimation of the OGWSBI algorithm and IOGWSBI;
FIG. 5 is a graph of root mean square error estimation and signal-to-noise ratio of the OGWSBI algorithm and the IOGWSBI algorithm DOA;
FIG. 6 is a graph of the relationship between the root mean square error of the OGWSBI algorithm and the DOA estimation of the IOGWSBI algorithm and the snapshot number;
FIG. 7 is a graph of the running time versus signal-to-noise ratio of the OGWSBI algorithm and the IOGWSBI algorithm.
Detailed Description
The technical solutions in the embodiments of the present invention will be clearly and completely described below with reference to the drawings in the embodiments of the present invention, and it is obvious that the described embodiments are only a part of the embodiments of the present invention, and not all of the embodiments. All other embodiments, which can be derived by a person skilled in the art from the embodiments given herein without making any creative effort, shall fall within the protection scope of the present invention.
As shown in fig. 1, a coherent source-exit-grating DOA estimation method based on weighted second-order sparse bayes is used for the following models: suppose K narrow-band far-field signals are at an angle [ theta ]1,θ2,...θK]TThe signal is incident on a uniform linear array consisting of M omnidirectional receiving array elements. The wavelength of the incident signal is lambda, and the distance between adjacent array elements is d which is lambda/2.
Figure BDA0003092317380000052
Is distributed in [0, π]And (3) discrete sampling values of the space angle, wherein N is the number of grid divisions. When an off-grid condition occurs, θ is due to the angle of the incident signalkNot falling completely on the grid of the observation matrix and therefore having an offset β from the angle on the observation matrix. Receive a signal of
Figure BDA0003092317380000061
The sparse model expression is YcPhi (beta) S + N, wherein
Figure BDA0003092317380000062
Is a signal source and is a signal source,
Figure BDA0003092317380000063
representing noise, wherein M is an array element number, N is a grid division number, T is a fast beat number, K is a signal source number, and phi is an over-complete dictionary formed by guide vectors; in the prior art, a first-order Taylor series approximation is adopted for a general guide vector;
the DOA estimation method comprises the following steps:
step S1: for received signal YcUnitary transformation real number domain conversion and signal de-coherent processing are carried out, so that a receiving signal Y after unitary transformation is obtained;
the signal Y is received in step S1cThe offset β at this time when unitary transformation is performed is 0, and the specific unitary transformation process is as follows:
step S11: for received signal YcConstructing an augmentation matrix to form a central Hermite matrix:
Figure BDA0003092317380000064
wherein ,Yc *Is YcThe complex conjugate of (a) and (b),
Figure BDA0003092317380000065
step S12: centering Hermite matrix YaugUnitary transformation processing:
Figure BDA0003092317380000066
wherein Y represents YaugAfter unitary transformation, transforming the complex number domain into a matrix of a real number domain, wherein U is a unitary matrix;
when M is even, UMCan be expressed as:
Figure BDA0003092317380000067
wherein, I is an M/2 xM/2 unit matrix, J is an M/2 xM/2 exchange matrix, the minor diagonal of the exchange matrix is 1, and other elements are 0;
when M is odd, UMCan be expressed as:
Figure BDA0003092317380000068
wherein, the dimension of I and J is (M-1)/2, I is an identity matrix of M/2 xM/2, J is a switching matrix of M/2 xM/2, the minor diagonal is 1, and other elements are 0.
In the same way, U2TIs a 2T-dimensional unitary matrix.
The unitary transformation can also recover the rank of the covariance matrix of the received signals, so that the rank is equal to the number K of the signal sources, i.e., rank (r) ═ K, and the decoherence performance is further realized.
Step S2: decomposing the singular value of Y of the received signal after unitary transformation to obtain YsvAnd SsvCombining the MUSIC algorithm to construct a weighting vector w; (ii) a
The method comprises the following specific steps: performing singular value decomposition on the received signal Y after unitary transformation, namely Y ═ Us∑sVs T+Ue∑eVe TTo obtain a signal subspace UsNoise subspace Ue(ii) a Order to
Figure BDA0003092317380000071
Figure BDA0003092317380000072
Figure BDA0003092317380000073
For the real signal source, Y is obtainedsv、Ssv；
The inverse MUSIC spatial spectrum formula is:
Figure BDA0003092317380000074
calculate the weight value w ═ w1,…,wN]TMeanwhile, normalizing the weight to obtain a weighted vector w: w ═ w/min (w).
Step S3: solving the hyperparameter alpha by combining the weighted sparse Bayesian inference to the unitary transformation of the parameter matrix of the overcomplete dictionary constructed by the second-order guide vector0Alpha, and offset beta, and outputs a DOA estimate
Figure BDA0003092317380000075
Step S3 includes the following steps:
step S31: constructing an over-complete dictionary by using a real source arrival angle guide vector of second-order Taylor expansion:
Figure BDA0003092317380000076
wherein ,AcFor the observation matrix:
Figure BDA0003092317380000077
Figure BDA0003092317380000078
is distributed in [0, π]Discrete sampling values of the space angles are obtained, and N is the number of grid divisions; the guide vector of the second-order Taylor expansion is
Figure BDA0003092317380000079
Figure BDA00030923173800000710
At an angle theta to the true signalkThe closest value of the grid to the grid value,
Figure BDA00030923173800000711
a (theta) is in
Figure BDA00030923173800000712
The first derivative of (a) is,
Figure BDA00030923173800000713
a (theta) is in
Figure BDA00030923173800000714
Second derivative of (a) of
Figure BDA00030923173800000715
The matrix is composed of
Figure BDA00030923173800000716
Figure BDA00030923173800000717
Composed matrix
Figure BDA00030923173800000718
Step S32: for A in overcomplete dictionaryc、BcAnd CcOver-complete dictionary after unitary transformation is obtained by performing unitary transformation
Figure BDA00030923173800000719
For uniform linear arrays, the pilot vectors satisfy the characteristics of center Hermite, and the step S32 is to compare A in overcomplete dictionaryc、BcAnd CcA, B and C for unitary transformation are shown below:
Figure BDA0003092317380000081
Figure BDA0003092317380000082
Figure BDA0003092317380000083
unitary transformation is carried out on the received signals and the observation matrix, so that the whole estimation model is converted from a complex number domain to a real number domain, and the calculated amount is effectively reduced.
Step S33: by using
Figure BDA0003092317380000084
Λ＝diag(α)、p(αn)＝Γ(α|1,wnRho) which is a user selectivity parameter and a signal source SsvA posterior probability distribution of
Figure BDA0003092317380000085
Figure BDA0003092317380000086
Obtaining the mean value and the variance mu (t) ═ alpha of the signal0∑ΦHysv(t),t＝1,....K，∑＝(α0ΦHΦ+Λ-1)-1Updating the mean value mu (t) and the variance sigma of the signal using the overcomplete dictionary phi after the unitary transformation, and updating the hyperparameter alpha based on the mean value mu (t) and the variance sigma of the signal0α and offset β; judging whether the maximum iteration number i is reachedmax or ||αi+1-αi||2/||αi||2If the two are not satisfied, updating the overcomplete dictionary phi after unitary transformation by using the updated offset beta, and continuously updating the mean value mu (t) and the variance sigma of the signal by using the updated overcomplete dictionary phi so as to update the hyper-parameter alpha0α and offset β; if one or two conditions are met, the iteration is ended, and the DOA estimated value is output
Figure BDA0003092317380000087
In step S33, the weighting vector w obtained in step S2 needs to be fused into a sparse bayesian framework. For signal distribution according to the sparse bayesian algorithm, it is assumed that it satisfies:
Figure BDA0003092317380000088
wherein ,ssv(t) represents SsvColumn t, diagonal matrix Λ ═ diag (α), α ═ α1,...,αN]T，αnControlling sparsity of signal source, i.e. only alpha in direction of presence of sourcenIs a non-zero number. The weight vector is added to its probability density hypothesis, assuming it satisfies the Gamma distribution: p (a)n)＝Γ(a|1,wnρ). Wherein Γ (·) is a Gamma function, ρ is a user-selective parameter, w is a user-selective parameternIs the nth weight coefficient. From a signal source SsvA posterior probability distribution of
Figure BDA0003092317380000089
The mean and covariance of the obtained signals satisfy:
μ(t)＝α0∑ΦHysv(t),t＝1,....K (2)
∑＝(α0ΦHΦ+Λ-1)-1 (3)
updated hyper-parameter alpha0And alpha is respectively:
Figure BDA0003092317380000091
Figure BDA0003092317380000092
in formula (4) and formula (5), μ ═ μ (1),.. μ (K)]＝α0∑ΦHYsv，
Figure BDA0003092317380000093
Figure BDA0003092317380000094
ρ＝ρ/K，
Figure BDA0003092317380000095
b,c→0。
The method for solving the angular deviation beta is to use the EM algorithm to make E { logp (Y)sv|Ssv,α0β) p (β) } max, minimizing the following:
Figure BDA0003092317380000096
the updated formula for the offset β is calculated as:
Figure BDA0003092317380000097
wherein ,
Figure BDA0003092317380000098
Figure BDA0003092317380000099
Figure BDA00030923173800000910
to further simplify the offset operation, let
Figure BDA00030923173800000911
In order to transform the matrix, the matrix is,
Figure BDA00030923173800000912
is the N × 1 column vector for row i 1 and the other rows 0. Beta is a-kRepresenting the part of the offset beta other than the kth element, and
Figure BDA00030923173800000913
for beta isnewTransforming to obtain:
Figure BDA0003092317380000101
wherein ,
Figure BDA0003092317380000102
d1 is of the formula betaTP1Beta is neutral to betakAn extraneous portion. Similarly:
V1 Tβ＝(V1 T)kβk+D2
Figure BDA0003092317380000103
Figure BDA0003092317380000104
Figure BDA0003092317380000105
wherein D2, D3, D4 and D5 are all betakThe extraneous part, then the offset βkIs the real root of the formula:
Figure BDA0003092317380000106
in addition, since the offset β has the same sparse structure as the source, P and V can be reduced to K or K × K dimensions. F' (β) in equation (7)k) The resulting root is the offset in one iteration, and during each iteration, the mean, covariance matrix, and the updated hyper-parameter α are calculated0α, β, and terminating the iteration condition to reach a maximum number of iterations imax or ||αi+1-αi||2/||αi||2＜τ。
Hyperparameter alpha0The updating process of α and the offset β is as follows:
(1) inputting an initial value: β is 0, ρ is 0.01, and c is d is 1 × 10-4,
Figure BDA0003092317380000107
Figure BDA0003092317380000108
(2) Obtaining an overcomplete dictionary phi after unitary transformation according to a formula (1), and obtaining a mean value mu (t) and a covariance sigma of a signal source according to formulas (2) and (3);
(3) obtaining the updated hyper-parameter alpha by combining the mean value mu (t) and the covariance sigma with the formula (4) and the formula (5)0α, updating the offset β using equation (7);
(4) judging whether the maximum iteration number i is reachedmax or ||αi+1-αi||2/||αi||2If both are not satisfied, updating the overcomplete dictionary phi after unitary transformation by using the updated offset beta according to the formula (1), namely turning to the step (1); if one or two conditions are met, the iteration is ended, and the DOA estimated value is output
Figure BDA0003092317380000111
Figure BDA0003092317380000112
wherein ,
Figure BDA0003092317380000113
denotes a mesh value, β 'of the j-th signal source'jThe grid offset for the jth signal source.
Next, a specific example is described, assuming that K is 2 narrow-band far-field signals at an angle θ1,θ2,...θK]TIncident on a uniform linear array consisting of 8 omnidirectional receiving array elements. The wavelength of the incident signal is lambda, the distance between adjacent array elements is d which is lambda/2,
Figure BDA0003092317380000114
is distributed in [0, π]The above space angle discrete sampling value, N is the number of grid divisions, N is 91, and the iteration number imax＝2000, margin parameter τ 10-3The specific parameters are shown in table 1:
TABLE 1 System simulation parameters
Parameter name Value of parameter
Array element number (M) 8
Array element interval (d) Lambda (wavelength)/2
Angular range 0°～180°
Angular interval (r) 2°
Information source number (K) 2
Fig. 2 is an angle section view of the OGWSBI algorithm at coherent and incoherent sources, which can be obtained from fig. 2, and the OGWSBI algorithm can accurately estimate the incoming direction of a signal when the signal source is an incoherent source; when the signal source is a coherent source, the estimation performance degrades rapidly.
Figure 3 is a cross-sectional view of the MUSIC algorithm estimating the angle of the coherent source before and after a unitary transform. As can be seen from fig. 3, after the unitary transform process, the rank of the covariance is restored to rank (r) ═ K, the signal subspace and the noise subspace are correctly represented, and the MUSIC algorithm can estimate the coherent source with high accuracy, so that the correct weighting vector can be provided.
The OGWSBI algorithm is improved based on the unitary transformation of the above second-order taylor expansion off-grid model of the guide vector and the estimation model, fig. 4 shows a comparison graph of the OGWSBI and the space spectrum of the monte carnot experiment of the algorithm IOGWSBI proposed herein at the next time when SNR is 10dB, and table 2 shows the corresponding results.
From fig. 4 and table 2, it can be seen that the estimated performance of IOGWSBI is significantly better than the OGWSBI algorithm. This shows that, when the signal source is coherent, the IOGWSBI algorithm can still maintain superior DOA estimation performance due to unitary transform decorrelation processing.
TABLE 2 DOA Settlement results
θ1 θ2
true signal 64.8° 86.5°
OGWSBI 65.06° 87°
IOGWSBI 64.88° 86.55°
Fig. 5 is a graph of the root mean square error of the DOA estimate as a function of the signal-to-noise ratio at a snapshot count L of 100. Fig. 6 is a plot of root mean square error of DOA estimates as a function of snapshot number at SNR of 10 dB. As can be seen from fig. 5 and 6, the root mean square error of both algorithms decreases gradually as the signal-to-noise ratio and the number of fast beats increase. But the weighting vector is correctly represented due to the decorrelation process performed by the IOGWSBI algorithm. Therefore, under the condition of equal signal-to-noise ratio or fast beat number, the proposed algorithm is superior to the OGWSBI algorithm, and the estimated performance is improved by about 50%.
The IOGWSBI algorithm is composed of 3 parts of unitary transformation, singular value decomposition and weighted sparse reconstruction. The unitary transform calculation is negligible since it involves only addition and subtraction operations. The computation complexity of the singular value decomposition is O (max (M2T, MT 2)). It can be found that the calculation only needs to be carried out once whether the unitary transformation or the singular value decomposition is carried out, while the algorithm is an iteration algorithm, and the complexity of the algorithm depends on the weighted sparse reconstruction part. The computational complexity of a single iteration of the OGWSBI algorithm is O (max (MN2)), and under a coherent source, the number of iterations increases due to wrong weighting information, further increasing the computational complexity. IOGWSBI uses unitary transformation, one real number field operation is 1/4 of previous complex number field, and the iteration number is greatly reduced due to correct weighting information. In addition, the IOGWSBI algorithm adopts a guide vector second-order Taylor expansion model, so that the operation complexity is improved to a certain extent compared with that of a first-order model, and the overall operation speed is still better than that of the OGWSBI algorithm.
Fig. 7 is a graph of the algorithm running time as a function of the signal-to-noise ratio at a snapshot count L of 100. As can be seen from FIG. 7, when the signal-to-noise ratio varies from 0dB to 20dB, the operation time of the two algorithms is continuously reduced along with the improvement of the signal-to-noise ratio, and the running efficiency of the IOGWSBI algorithm is improved by about 60% compared with that of the OGWSBI algorithm. Therefore, under the condition of a coherent information source, the IOGWSBI algorithm is an efficient and high-precision off-grid DOA estimation algorithm and is more suitable for scenes with higher precision requirements.
On the basis of the OGWSBI algorithm, unitary transformation processing is adopted for the estimation model, so that coherent information sources can be subjected to coherent elimination processing, the problem of insufficient estimation performance of the OGWSBI algorithm under the coherent information sources is solved, the algorithm can be converted from a complex number domain to a real number domain for operation, and the operation efficiency of the algorithm is remarkably improved. In addition, in order to further improve the estimation accuracy, second-order taylor expansion is also performed on the basis of first-order taylor expansion of the guide vector, so that the fitting error of the signal is reduced. The algorithm provided by the invention is effectively improved in precision and efficiency, and the actual application scene of the weighted sparse Bayesian algorithm is widened.
The above description is only an embodiment of the present invention, and is not limited thereto, and any changes or substitutions that are within the spirit and principle of the present invention are included in the protection scope of the present invention, and therefore, the protection scope of the present invention should be subject to the protection scope of the claims.

Claims (7)

1. A coherent source-off-grid DOA estimation method based on weighted second-order sparse Bayes is characterized in that the DOA estimation method is used for the following models: receive a signal of
Figure FDA0003092317370000011
Its sparse expression is YcPhi (beta) S + N, wherein
Figure FDA0003092317370000012
Is a signal source, beta is an angular deviation,
Figure FDA0003092317370000013
representing noise, wherein M is an array element number, N is a grid division number, T is a fast beat number, K is a signal source number, and phi is an over-complete dictionary formed by guide vectors;
the DOA estimation method comprises the following steps:
step S1: for received signal YcUnitary transformation realizes real number domain transformation and signal decorrelation processing, thereby obtaining a receiving signal Y after unitary transformation;
step S2: decomposing the singular value of Y of the received signal after unitary transformation to obtain YsvAnd SsvCombining the MUSIC algorithm to construct a weighting vector w;
step S3: solving the hyperparameter alpha by combining the weighted sparse Bayesian inference to the unitary transformation of the parameter matrix of the overcomplete dictionary constructed by the second-order guide vector0Alpha, and offset beta, and outputs a DOA estimate
Figure FDA0003092317370000014
2. The method as claimed in claim 1, wherein the step S3 comprises the following steps:
step S31: constructing an over-complete dictionary by using a real source arrival angle guide vector of second-order Taylor expansion:
Figure FDA0003092317370000015
wherein ,AcFor the observation matrix:
Figure FDA0003092317370000016
is distributed in [0, π]Discrete sampling values of the space angles are obtained, and N is the number of grid divisions; the guide vector of the second-order Taylor expansion is
Figure FDA0003092317370000017
Figure FDA0003092317370000018
At an angle theta to the true signalkThe closest value of the grid to the grid value,
Figure FDA0003092317370000019
a (theta) is in
Figure FDA00030923173700000110
The first derivative of (a) is,
Figure FDA00030923173700000111
a (theta) is in
Figure FDA00030923173700000112
Second derivative of (a) of
Figure FDA00030923173700000113
The matrix is composed of
Figure FDA00030923173700000114
Figure FDA00030923173700000115
Composed matrix
Figure FDA00030923173700000116
Step S32: for A in overcomplete dictionaryc、BcAnd CcOver-complete dictionary after unitary transformation is obtained by performing unitary transformation
Figure FDA00030923173700000117
Step S33: by using
Figure FDA0003092317370000021
Λ＝diag(α)、p(αn)＝Γ(α|1,wnRho) which is a user selectivity parameter and a signal source SsvA posterior probability distribution of
Figure FDA0003092317370000022
Figure FDA0003092317370000023
Obtaining the mean value and the variance mu (t) ═ alpha of the signal0∑ΦHysv(t),t＝1,....K，∑＝(α0ΦHΦ+Λ-1)-1Updating the mean value mu (t) and the variance sigma of the signal using the overcomplete dictionary phi after the unitary transformation, and updating the hyperparameter alpha based on the mean value mu (t) and the variance sigma of the signal0Alpha and offsetBeta; judging whether the maximum iteration number i is reachedmax or ||αi+1-αi||2/||αi||2If the two are not satisfied, updating the overcomplete dictionary phi by using the updated offset beta, and continuously updating the mean value mu (t) and the variance sigma of the signal by using the updated overcomplete dictionary phi so as to update the overcompensate parameter alpha0α and offset β; if one or two conditions are met, the iteration is ended, and the DOA estimated value is output
Figure FDA0003092317370000024
3. The method as claimed in claim 2, wherein the updated hyper-parameter α is0And alpha is respectively:
Figure FDA0003092317370000025
Figure FDA0003092317370000026
wherein, μ ═ μ (1),.. μ (K)]＝α0∑ΦHYsv，
Figure FDA0003092317370000027
ρ＝ρ/K，
Figure FDA0003092317370000028
b,c→0。
4. The method as claimed in claim 3, wherein the updated offset β is:
Figure FDA0003092317370000029
wherein ,
Figure FDA00030923173700000210
Figure FDA00030923173700000211
Figure FDA0003092317370000031
5. the method as claimed in claim 2, wherein the step S1 comprises the following steps:
step S11: for received signal YcConstructing an augmentation matrix to form a central Hermite matrix:
Figure FDA0003092317370000039
wherein ,Yc *Is YcThe complex conjugate of (a) and (b),
Figure FDA0003092317370000032
step S12: centering Hermite matrix YaugUnitary transformation processing:
Figure FDA0003092317370000033
wherein Y represents YaugAfter unitary transformation, transforming the complex number domain into a matrix of a real number domain, wherein U is a unitary matrix;
when M is even, UMCan be expressed as:
Figure FDA0003092317370000034
wherein, I is an M/2 xM/2 unit matrix, J is an M/2 xM/2 exchange matrix, the minor diagonal of the exchange matrix is 1, and other elements are 0;
when M is odd, UMCan be expressed as:
Figure FDA0003092317370000035
wherein the dimensions of I and J are (M-1)/2;
in the same way, U2TIs a 2T-dimensional unitary matrix.
6. The method as claimed in claim 5, wherein the step S32 is performed on A in the overcomplete dictionary based on the weighted second-order sparse Bayes coherent source-exit-grating DOA estimation methodc、BcAnd CcA, B and C for unitary transformation are shown below:
Figure FDA0003092317370000036
Figure FDA0003092317370000037
Figure FDA0003092317370000038
7. the method as claimed in claim 2, wherein the step S2 comprises the following steps:
performing singular value decomposition on the received signal Y after unitary transformation, namely Y ═ Us∑sVs T+Ue∑eVe TTo obtain a signal subspace UsNoise subspace Ue(ii) a Order to
Figure FDA0003092317370000041
Figure FDA0003092317370000042
For the real signal source, Y is obtainedsv、Ssv；
The inverse MUSIC spatial spectrum formula is:
Figure FDA0003092317370000043
calculate the weight value w ═ w1,…,wN]TMeanwhile, normalizing the weight to obtain a weighted vector w: w ═ w/min (w).
CN202110599238.0A 2021-05-31 2021-05-31 Coherent source grid-off DOA estimation method based on weighted second-order sparse Bayes Active CN113534040B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202110599238.0A CN113534040B (en) 2021-05-31 2021-05-31 Coherent source grid-off DOA estimation method based on weighted second-order sparse Bayes

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202110599238.0A CN113534040B (en) 2021-05-31 2021-05-31 Coherent source grid-off DOA estimation method based on weighted second-order sparse Bayes

Publications (2)

Publication Number Publication Date
CN113534040A true CN113534040A (en) 2021-10-22
CN113534040B CN113534040B (en) 2023-08-11

Family

ID=78094898

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202110599238.0A Active CN113534040B (en) 2021-05-31 2021-05-31 Coherent source grid-off DOA estimation method based on weighted second-order sparse Bayes

Country Status (1)

Country Link
CN (1) CN113534040B (en)

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116626646A (en) * 2023-07-21 2023-08-22 西安电子科技大学 Radar target gridding-free loss coherent accumulation method based on time-frequency non-uniform sampling

Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
KR20080086950A (en) * 2007-03-23 2008-09-29 삼성전자주식회사 Method and apparatus for coherent source doa estimation
CN109116293A (en) * 2018-08-22 2019-01-01 上海师范大学 A kind of Wave arrival direction estimating method based on sparse Bayesian out of place
CN109490819A (en) * 2018-11-16 2019-03-19 南京邮电大学 A kind of Wave arrival direction estimating method out of place based on management loading
CN110109050A (en) * 2019-04-08 2019-08-09 电子科技大学 The DOA estimation method of unknown mutual coupling under nested array based on sparse Bayesian
CN111337893A (en) * 2019-12-19 2020-06-26 江苏大学 Off-grid DOA estimation method based on real-value sparse Bayesian learning

Patent Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
KR20080086950A (en) * 2007-03-23 2008-09-29 삼성전자주식회사 Method and apparatus for coherent source doa estimation
CN109116293A (en) * 2018-08-22 2019-01-01 上海师范大学 A kind of Wave arrival direction estimating method based on sparse Bayesian out of place
CN109490819A (en) * 2018-11-16 2019-03-19 南京邮电大学 A kind of Wave arrival direction estimating method out of place based on management loading
CN110109050A (en) * 2019-04-08 2019-08-09 电子科技大学 The DOA estimation method of unknown mutual coupling under nested array based on sparse Bayesian
CN111337893A (en) * 2019-12-19 2020-06-26 江苏大学 Off-grid DOA estimation method based on real-value sparse Bayesian learning

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
J YANG ET AL.: "An efficient off-grid DOA estimation approach for nested array signal prosessing by using sparse Bayesian learning strategies", SIGNAL PROCESSING *
高阳;陈俊丽;杨广立;: "基于酉变换和稀疏贝叶斯学习的离格DOA估计", 通信学报, no. 06 *

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116626646A (en) * 2023-07-21 2023-08-22 西安电子科技大学 Radar target gridding-free loss coherent accumulation method based on time-frequency non-uniform sampling
CN116626646B (en) * 2023-07-21 2023-09-22 西安电子科技大学 Radar target gridding-free loss coherent accumulation method based on time-frequency non-uniform sampling

Also Published As

Publication number Publication date
CN113534040B (en) 2023-08-11

Similar Documents

Publication Publication Date Title
CN109444810B (en) Mutual-prime array non-grid DOA estimation method under nonnegative sparse Bayesian learning framework
CN110208735B (en) Sparse Bayesian learning-based coherent signal DOA estimation method
CN103971029B (en) Alternant iteration method for DOA (direction of arrival) estimation under grid mismatch
Liu et al. An efficient maximum likelihood method for direction-of-arrival estimation via sparse Bayesian learning
CN110109050B (en) Unknown mutual coupling DOA estimation method based on sparse Bayes under nested array
CN109116293B (en) Direction-of-arrival estimation method based on lattice-separated sparse Bayes
CN111337893A (en) Off-grid DOA estimation method based on real-value sparse Bayesian learning
CN110244272B (en) Direction-of-arrival estimation method based on rank-denoising model
CN110749855B (en) Covariance domain nulling-based uniform linear array direction-of-arrival estimation method
CN113534040B (en) Coherent source grid-off DOA estimation method based on weighted second-order sparse Bayes
Trinh-Hoang et al. An improved DOA estimator based on partial relaxation approach
CN113376569B (en) Nested array sparse representation direction-of-arrival estimation method based on maximum likelihood
CN112731273B (en) Low-complexity signal direction-of-arrival estimation method based on sparse Bayesian
CN113567913A (en) Two-dimensional plane DOA estimation method based on iteration reweighting dimension reduction
CN110174657B (en) Direction-of-arrival estimation method based on rank-one dimension reduction model and block matrix recovery
CN112924947A (en) MIMO radar steady DOA estimation method based on real-value sparse Bayesian learning
CN116299193A (en) MIMO radar intelligent DOA estimation method
CN115587281A (en) Array element failure MIMO radar angle estimation method based on factor matrix prior
CN115130504A (en) Robust beam forming method based on sparse Bayesian learning
CN113406560A (en) Angle and frequency parameter estimation method of incoherent distributed broadband source
CN111077493B (en) Neural array direction of arrival estimation method based on real-value lattice separation variation Bayesian inference
Li et al. DOA Estimation Based on Sparse Reconstruction via Acoustic Vector Sensor Array under Non-uniform Noise
CN114624665B (en) Mutual coupling error DOA self-correction method based on dynamic parameter iterative optimization
CN115422732B (en) Mutual coupling optimization array, design method thereof and coherent signal AOA estimation method
CN114415110B (en) Direct positioning method for non-negative sparse Bayesian learning

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant