CN112731273B - A low-complexity signal direction-of-arrival estimation method based on sparse Bayesian - Google Patents

Abstract

The invention discloses a sparse Bayesian-based low-complexity signal direction-of-arrival estimation method, and aims to solve the technical problems of high computational complexity and low computational efficiency of an MMV model in the prior art when time correlation problems are considered. It comprises the following steps: expanding an incident signal, and obtaining the prior of the expanded incident signal by using a Markov probability prior model; and (3) obtaining an input function and an output function of a GAMP algorithm according to the priori of the expanded incident signal, performing GAMP iteration, obtaining a recovered expanded incident signal, and obtaining a signal direction of arrival according to the recovered expanded incident signal. The method can reduce the calculation complexity of DOA estimation on the premise of considering time-related factors, improve the calculation efficiency and accurately estimate the signal arrival direction.

Description

A low-complexity signal direction-of-arrival estimation method based on sparse Bayes

Technical Field

The present invention relates to a low-complexity signal direction-of-arrival estimation method based on sparse Bayes, and belongs to the technical field of signal processing.

Background Art

Direction of Arrival (DOA) is the estimation of the direction of arrival of an incident signal using the received data of an antenna array in a noisy environment. Its basic principle is to use the phase difference between the received data of different array elements in the spatial array to estimate the spatial angle at which the signal arrives at the array. It can be applied to radio communications, radar, sonar, navigation, earthquake detection, biomedicine and other fields, and is of great significance.

The SBL algorithm, which has emerged in recent years, was originally proposed by Tipping et al. as a machine learning algorithm around 2001, and was subsequently introduced into the field of sparse signal recovery/compressed sensing. At first, SBL was only applied to the model of a single measurement vector (SMV), and later gradually expanded to multiple measurement vectors (MMV). It has the advantages of high resolution and reduced computational complexity, but the DOA model used in MMV is relatively idealized and does not consider more practical factors. In 2011, Zhang Zhilin et al. applied SBL to the time-related direction finding scenario and proposed the Temporally Sparse Bayesian Learning (TSBL) algorithm, introducing a hyperparameter to control the time correlation between snapshots. Although this can solve the time correlation problem, the introduced hyperparameter will increase the computational complexity and make the efficiency very low. If there are large-scale problems such as a large number of array elements, accurate direction finding cannot be achieved, and further algorithm optimization is required.

Summary of the invention

In order to solve the problems of high computational complexity and low computational efficiency of the MMV model in the prior art when considering time-related issues, the present invention proposes a low-complexity signal direction of arrival estimation method based on sparse Bayes. On the basis of considering time-related factors, the Markov prior probability distribution is used to optimize the derivation process of obtaining the posterior probability, and the GAMP algorithm is used for decoupling to replace the matrix inversion step of the prior art, thereby reducing the amount of calculation and improving the computational efficiency.

In order to solve the above technical problems, the present invention adopts the following technical means:

The present invention proposes a low-complexity signal direction of arrival estimation method based on sparse Bayes, comprising the following steps:

Expanding the incident signal based on spatial domain grid partitioning;

The Markov probability prior model is used to obtain the prior of the expanded incident signal;

The input function and output function of the GAMP algorithm are obtained according to the prior of the expanded incident signal;

Using the input function and output function of the GAMP algorithm, the GAMP iteration is performed on the snapshot of the expanded incident signal to obtain the restored expanded incident signal;

The signal arrival direction is obtained according to the restored and expanded incident signal.

Furthermore, the specific steps of incident signal expansion are as follows:

Assume that the incident signal is x, and the receiving array is composed of M uniform linear array antennas. When the receiving array receives T snapshots of the incident signal, the array output model is obtained:

y＝A(θ)x+e (1)

Where y is the received signal of the receiving array, A(θ) is the array flow matrix, and e is the Gaussian white noise received by the receiving array;

其中，

表示第n个超完备的角度，H为超完 备的角度集合

中的角度数量，n＝1,2,…,H；The spatial domain of the receiving array is gridded using the exhaustive method to obtain an overcomplete angle set.

in,

represents the nth supercomplete angle, H is the supercomplete angle set

The number of angles in, n = 1, 2, ..., H;

和A(θ)获得扩展后的阵列流型矩阵：based on

and A(θ) to obtain the expanded array flow matrix:

Among them, A(θ)_sparse represents the expanded array flow matrix;

The overcomplete array output model is obtained based on the expanded array flow matrix:

表示扩展后的接收信号，

表示扩展后的入射信号。in,

represents the received signal after expansion,

represents the incident signal after expansion.

Furthermore, the spatial domain of the receiving array is [-90°, 90°].

Furthermore, the steps for obtaining the prior of the expanded incident signal are as follows:

The expanded incident signal is modeled using the Markov probability prior model to obtain the relationship between the nth expanded incident signal at the tth snapshot and the t-1th snapshot:

表示第t个快拍的第n个扩展后的入射信号，

表示第t-1 个快拍的第n个扩展后的入射信号，

表示

和

之间的关系， β为时间相关系数，β∈(-1,1)，γ_n为第n个扩展后的入射信号的先验误差， t＝1,2,…,T；in,

represents the nth expanded incident signal of the tth snapshot,

represents the nth expanded incident signal of the t-1th snapshot,

express

and

, β is the time correlation coefficient, β∈(-1,1), γ _n is the prior error of the incident signal after the nth expansion, t＝1,2,…,T;

和后一快拍的影响

其中，

表示入射信号

的先验概率，

表示第n个扩展后的入射信号t快 拍时t-1快拍前向传递的均值，

表示第n个扩展后的入射信号t快拍时 t-1快拍前向传递的方差，

表示入射信号

的先验概率，

表示第 n个扩展后的入射信号t快拍时t+1快拍后向传递的均值，

表示第n个扩 展后的入射信号t快拍时t+1快拍后向传递的方差；According to the relationship between the previous and next snapshots of the expanded incident signal, the influence of the previous snapshot is obtained.

And the impact of the next snapshot

in,

Represents the incident signal

The prior probability of

represents the mean of the forward pass of snapshot t-1 at the time of snapshot t of the nth extended incident signal,

represents the variance of the forward transmission of the t-1 snapshot at the t snapshot of the n-th expanded incident signal,

Represents the incident signal

The prior probability of

represents the mean value of the backward transmission of the t+1 snapshot of the incident signal after the n-th expansion,

represents the variance of the incident signal t+1 snapshot transmitted backward at the time of snapshot t after the nth expansion;

和

获得扩展后的 入射信号

的先验：according to

and

Get the expanded incident signal

Priors:

的表达式如 下：Furthermore, the nth expanded incident signal of the tth snapshot

The expression is as follows:

in,

Furthermore, the expression of the input function of the GAMP algorithm is as follows:

表示q^(t)的噪声方差，y^(t)表示第t个快拍接收矩阵的接收信号， σ²表示y^(t)的噪声方差；Where _gs represents the input function, q ^(t) represents the approximate value of the received signal of the t-th snapshot without noise influence,

represents the noise variance of q ^(t) , y ^(t) represents the received signal of the t-th snapshot receiving matrix, ^σ2 represents the noise variance of y ^(t) ;

The expression of the output function of the GAMP algorithm is as follows:

表示第t个快拍的第n个扩展后的入射信 号的近似值，

表示

的噪声方差。Among them, g _x represents the output function,

represents the approximate value of the n-th expanded incident signal of the t-th snapshot,

express

The noise variance.

Furthermore, the specific operation of obtaining the restored expanded incident signal is as follows:

Initialize the expanded incident signal of the t-th snapshot

进行GAMP迭代，更新 扩展后的入射信号的均值和方差；Using the input function and output function of the GAMP algorithm

Perform GAMP iteration to update the mean and variance of the expanded incident signal;

In each iteration process, the EM algorithm is used to update the iteration parameters, and the iteration convergence condition is used to judge the iteration convergence;

Obtaining a restored extended incident signal according to the mean and variance of the extended incident signal that satisfy the iterative convergence judgment;

Wherein, the iterative convergence condition includes the GAMP iterative convergence condition and the EM iterative convergence condition.

Furthermore, assuming that the maximum number of iterations of the GAMP algorithm is G, the iterative process of each iteration is as follows:

的方差

均值X^i(t)、先验误差γⁱ和接收信号的噪声方差(σ²)ⁱ， 其中，i∈[1,G]；According to the output of the i-1th iteration, the extended incident signal of the tth snapshot in the i-th iteration is obtained.

Variance

Mean Xi ^(t) , prior error γ ⁱ and noise variance of received signal (σ ² ) ⁱ , where i∈[1,G];

和X^i(t)计算第i次迭代中无噪声影响的第t 个快拍的接收信号的近似值q^i(t)：Let S = |A(θ)_sparse| ² , and use

The approximation qi ^(t) of the received signal of the t-th snapshot without noise influence in the ith iteration is calculated by using Xi ^(t) :

Where g _s ^(i-1)(t) represents the input function of the t snapshot of the i-1th iteration;

Update the input function g _s ^i(t) of the t snapshot of the i-th iteration using γ ⁱ , (σ ² ) ⁱ and qi ^(t) ;

Use the updated input function g _si ^(t) to calculate the approximate value ri ^(t) of the extended incident signal and its noise variance of the t-th snapshot.

为输入阻尼系数，

s^(i-1)(t)表示第i-1次迭代的值，S^T表示S的转 置，(g_s ^i(t))'表示输入函数g_s ^i(t)的导数；Among them, A(θ)_sparse ^T represents the transpose of A(θ)_sparse,

is the input damping coefficient,

s ^(i-1)(t) represents the value of the i-1th iteration, S ^T represents the transpose of S, and (g _s ^i(t) )' represents the derivative of the input function g _s ^i(t) ;

更新第i次迭代的t快拍的输出函数g_x ^i(t)；Using γ ⁱ , ri ^(t) and

Update the output function _gxi ^(t) of the t snapshot of the i-th iteration;

Update the variance and mean of the expanded incident signal according to the updated output function _gxi ^(t) :

表示第i次迭代输出的扩展后的入射信号的方差，(g_x ^i(t))'表 示输出函数g_x ^i(t)的导数，X^(i+1)(t)表示第i次迭代输出的扩展后的入射信号的 均值，

为输出阻尼系数，

in,

represents the variance of the expanded incident signal output at the i-th iteration, (g _x ^i(t) )' represents the derivative of the output function g _x ^i(t) , X ^(i+1)(t) represents the mean of the expanded incident signal output at the i-th iteration,

is the output damping coefficient,

进行GAMP迭代收敛判断， 满足GAMP迭代收敛条件时，结束迭代，输出第i次迭代的扩展后的入射 信号的均值和方差，不满足GAMP迭代收敛条件时，进入下一步；Using the GAMP iterative convergence condition, we can calculate X ^(i+1)(t) and

Perform GAMP iteration convergence judgment. When the GAMP iteration convergence condition is met, the iteration is terminated and the mean and variance of the incident signal after expansion of the i-th iteration are output. When the GAMP iteration convergence condition is not met, proceed to the next step.

分别计算第i+1次 迭代中的先验误差γⁱ⁺¹和噪声方差(σ²)ⁱ⁺¹，计算公式如下：Based on the EM algorithm, using the output of the i-th iteration X ^(i+1)(t) and

The prior error γ ⁱ⁺¹ and the noise variance (σ ² ) ⁱ⁺¹ in the i+1th iteration are calculated respectively, and the calculation formulas are as follows:

(σ ² ) ⁱ⁺¹ =||y ^(t) -A(θ)_sparse·X ^(i+1)(t) || ² (15)

进行EM迭代收敛判断，满足 EM迭代收敛条件时，输出第i+1次迭代的扩展后的入射信号的均值和方差， 结束迭代，不满足EM迭代收敛条件时，继续迭代。Using the EM iteration convergence condition, we can calculate X ^(i+1)(t) and

Perform EM iteration convergence judgment. When the EM iteration convergence condition is met, output the mean and variance of the expanded incident signal of the i+1th iteration, and end the iteration. When the EM iteration convergence condition is not met, continue the iteration.

Furthermore, the GAMP iteration convergence condition is:

或迭代次数i达到GAMP算法的最 大迭代次数G，其中，ε_gamp表示GAMP算法迭代的归一化公差参数。

Or the number of iterations i reaches the maximum number of iterations G of the GAMP algorithm, where ε _gamp represents the normalized tolerance parameter of the GAMP algorithm iteration.

Furthermore, the EM iteration convergence condition is:

或迭代次数i达到EM算法的最大迭 代次数K，其中，ε_em表示EM算法迭代的归一化公差参数。

Or the number of iterations i reaches the maximum number of iterations K of the EM algorithm, where ε _em represents the normalized tolerance parameter of the EM algorithm iteration.

The following advantages can be obtained by using the above technical means:

The present invention proposes a low-complexity signal direction of arrival estimation method based on sparse Bayesian, which uses Markov probability prior model to build a model, takes into account the time-related factors of the direction of arrival estimation in multiple blocks, and improves the accuracy of the DOA model under the large-block and time-related MMV model. The method of the present invention is suitable for DOA model estimation in non-ideal environments, and improves the robustness of DOA estimation under low signal-to-noise ratio conditions. After obtaining the prior of the incident signal, the present invention combines the GAMP idea, introduces intermediate parameters to perform GAMP iteration, and replaces the matrix inversion step when calculating the posterior probability density under the traditional Bayesian framework through iterative operation, which greatly reduces the computational complexity, improves the computational efficiency, and improves the resolution of the signal direction of arrival estimation.

In the process of signal arrival direction estimation, the present invention not only takes into account non-ideal environments such as time correlation, but also balances good performance such as high resolution and low complexity, breaking through the limitations of the existing SBL algorithm in terms of accuracy, complexity, robustness, etc., and is of great significance to the research of modern application scenarios.

BRIEF DESCRIPTION OF THE DRAWINGS

FIG1 is a flow chart of the steps of a low-complexity signal direction-of-arrival estimation method based on sparse Bayesian in the present invention.

FIG. 2 is a schematic diagram of message transmission between snapshots during multi-block shooting according to an embodiment of the present invention.

FIG3 is a normalized power spectrum diagram of the method of the present invention and the TMSBL algorithm under low signal-to-noise ratio conditions in an embodiment of the present invention.

FIG4 is a normalized power spectrum diagram of the method of the present invention and the TMSBL algorithm under a high signal-to-noise ratio in an embodiment of the present invention.

FIG5 is a comparison diagram of the RMSE between the method of the present invention and the TMSBL algorithm in an embodiment of the present invention.

FIG6 is a comparison diagram of the method of the present invention and the TMSBL algorithm in CUP-time in an embodiment of the present invention.

DETAILED DESCRIPTION

The technical solution of the present invention is further described below in conjunction with the accompanying drawings:

The present invention proposes a low-complexity signal direction of arrival estimation method based on sparse Bayes, as shown in FIG1 , which specifically includes the following steps:

Step 1: Expand the incident signal based on spatial domain grid division;

Step 2: Using the Markov probability prior model to obtain the prior of the expanded incident signal;

Step 3, obtaining the input function and output function of the GAMP algorithm according to the prior of the expanded incident signal;

Step 4: perform GAMP iteration on the snapshot of the expanded incident signal using the input function and output function of the GAMP algorithm to obtain a restored expanded incident signal;

Step 5: Obtain the signal arrival direction according to the restored and expanded incident signal.

In the embodiment of the present invention, W narrowband far-field signal sources with a wavelength of λ are provided, which are incident on the receiving end at an angle of _θw , w∈{1,2,···,W}, the incident signal is x, and the receiving array of the receiving end is composed of M uniform linear array antennas.

When the receiving array receives T snapshots of the incident signal, the array output model is obtained:

y＝A(θ)x+e (16)

Wherein, y is the received signal of the receiving array, A(θ) is the array flow matrix, e is the Gaussian white noise received by the receiving array, e satisfies the normal distribution, e～N(0,σ ² I), σ ² represents the noise variance of the received signal.

The expression of the array flow matrix is as follows:

对入射信号进行 扩展，其中，

表示第n个超完备的角度，H为超完备的角度集合

中的角 度数量，n＝1,2,…,H。Based on the sparse representation theory, the direction of the incident signal is sparse in the entire spatial domain, and the resolvable spatial domain of the receiving array is [-90°, 90°]. The spatial domain of the receiving array is gridded using the exhaustive method to obtain an over-complete angle set.

The incident signal is expanded, where

represents the nth supercomplete angle, H is the supercomplete angle set

The number of angles in, n = 1, 2, …, H.

和A(θ)获得扩展后的阵列流型矩阵A(θ)_sparse：based on

and A(θ) to obtain the expanded array flow matrix A(θ)_sparse:

Because the incident signal falls on the grid, the overcomplete array output model is obtained according to the expanded array flow matrix:

表示扩展后的接收信号，

表示扩展后的入射信号，

in,

represents the received signal after expansion,

represents the incident signal after expansion,

再利 用扩展后的入射信号

得到入射信号x的DOA方向。The present invention uses the received signal and A(θ)_sparse to recover the expanded incident signal.

Reusing the expanded incident signal

Get the DOA direction of the incident signal x.

In step 2, the steps for obtaining the prior of the expanded incident signal are as follows:

Step 201: Under the Bayesian framework, taking into account the time-related factors, that is, the incident signal at the current moment will be affected by the previous and next snapshots, the signal source is modeled using the Markov probability prior model, and the product of the two messages of the previous and next snapshots is combined into one message.

可 以表示为：The present invention introduces a time correlation coefficient β, β∈(-1,1), which can represent the relationship between the previous snapshot and the next snapshot of a signal. The incident signal after the nth expansion of the tth snapshot

It can be expressed as:

表示第t-1个快拍的第n个扩展后的入射信号，

γ_n为第n个扩展后的入射信号的先验误差，t＝1,2,…,T。in,

represents the nth expanded incident signal of the t-1th snapshot,

_γn is the a priori error of the incident signal after the nth expansion, t=1,2,…,T.

According to formula (19), the relationship between the incident signal after the nth expansion at the tth snapshot and the t-1th snapshot is obtained:

表示

和

之间的关系。in,

express

and

The relationship between.

为第t个快拍的第n个扩展后的入射信号的先 验，则多块拍时快拍之间消息传递如图2所示。Under the MMV model, let

is the prior of the nth extended incident signal of the tth snapshot, then the message transmission between snapshots during multi-block shooting is shown in Figure 2.

设后向消息传递的先 验也是一个高斯分布，引入参数υ和φ作为其均值和方差，则后一快拍的影 响为

其中，

表示入射信号

的先验概率，

表示第n个扩展后的入射信号t快拍时t-1快拍前向传递的均值，

表示第 n个扩展后的入射信号t快拍时t-1快拍前向传递的方差，

表示入射 信号

的先验概率，

表示第n个扩展后的入射信号t快拍时t+1快拍后 向传递的均值，

表示第n个扩展后的入射信号t快拍时t+1快拍后向传 递的方差。Assume that the prior of the forward message transmission is a Gaussian distribution, and introduce parameters η and ψ as its mean and variance, then the impact of the previous snapshot is

Assume that the prior of the backward message transmission is also a Gaussian distribution, and introduce parameters υ and φ as its mean and variance, then the impact of the next snapshot is

in,

Represents the incident signal

The prior probability of

represents the mean of the forward pass of snapshot t-1 at the time of snapshot t of the nth extended incident signal,

represents the variance of the forward transmission of the t-1 snapshot at the t snapshot of the n-th expanded incident signal,

Represents the incident signal

The prior probability of

represents the mean value of the backward transmission of the t+1 snapshot of the incident signal after the n-th expansion,

It represents the variance of the incident signal after the n-th expansion when the t+1 snapshot is transmitted backward.

和

获得第t个快 拍的第n个扩展后的入射信号

的先验：merge

and

Get the nth extended incident signal of the tth snapshot

Priors:

After obtaining the prior of the incident signal, the present invention combines the GAMP concept and performs subsequent calculations through intermediate parameters, replacing the related steps of matrix inversion, thereby reducing the computational complexity. The present invention introduces an intermediate parameter r as an approximation of the extended incident signal, given the noise variance τ _r of r, and introduces an intermediate parameter q as an approximation of A(θ)x, that is, an approximation of the received signal without noise influence, given the noise variance τ _q of q.

和方差

的表达式：According to the Gaussian probability density function and convolution operation, the mean of the forward message passing prior is obtained by substituting formula (20) into

and variance

The expression is:

so:

表示第t-1个快拍的第n个扩展后的入射信号的近似值，

表示

的噪声方差。in,

represents the approximate value of the n-th expanded incident signal of the t-1-th snapshot,

express

The noise variance.

和方差

的表达式：Similarly, the mean of the backward message passing prior can be obtained

and variance

The expression is:

According to the definition of input function in GAMP algorithm, the expression of input function can be derived:

表示q^(t)的噪声方差，z＝A(θ)x，y^(t)表示第t个快拍接收矩阵 的接收信号，σ²表示y^(t)的噪声方差。Where _gs represents the input function, q ^(t) represents the approximate value of the received signal of the t-th snapshot without noise influence,

represents the noise variance of q ^(t) , z=A(θ)x, y ^(t) represents the received signal of the t-th snapshot receiving matrix, ^{and σ2} represents the noise variance of y ^(t) .

可以推导出输出函数 的表达式：According to the definition and prior of the output function in the GAMP algorithm

The expression of the output function can be derived:

表示第t个快拍的第n个扩展后的入射信 号的近似值，

表示

的噪声方差。Among them, g _x represents the output function,

represents the approximate value of the n-th expanded incident signal of the t-th snapshot,

express

The noise variance.

In the embodiment of the present invention, the specific operations of step 4 are as follows:

Initialize the expanded incident signal of the t-th snapshot

进行GAMP迭代，更新 扩展后的入射信号的均值和方差；Using the input function and output function of the GAMP algorithm

Perform GAMP iteration to update the mean and variance of the expanded incident signal;

In each iteration process, the EM algorithm is used to update the iteration parameters, and the iteration convergence condition is used to judge the iteration convergence;

Obtaining a restored extended incident signal according to the mean and variance of the extended incident signal that satisfy the iterative convergence judgment;

Among them, the iterative convergence conditions include the GAMP iterative convergence conditions and the EM iterative convergence conditions. In each iteration process, two convergence judgments must be made based on the two iterative convergence conditions. If any one of the convergence judgments is passed, the iteration will terminate.

Assume that the maximum number of iterations of the GAMP algorithm is G. When performing the first iteration, it is necessary to initialize the number of iterations, the mean, variance, prior error, and noise variance of the extended incident signal. As the number of iterations increases, the above values will also change until the iterative convergence conditions are met and the restored extended incident signal is output.

Taking the i-th iteration as an example, the iteration process is as follows:

的方差

均值X^i(t)、先验误差γⁱ和接收信号的噪声方差 (σ²)ⁱ，其中，i∈[1,G]。(1) Based on the output of the i-1th iteration, obtain the extended incident signal of the tth snapshot in the i-th iteration

Variance

mean Xi ^(t) , prior error γ ⁱ and noise variance of received signal (σ ² ) ⁱ , where i∈[1,G].

(2) Let S = |A(θ)_sparse| ² , and obtain the calculation formula of the variance of the approximate value of the received signal:

表示第i次迭代中无噪声影响的第t个快拍的接收信号的近 似值的方差。in,

represents the variance of the approximation of the received signal of the t-th snapshot in the ith iteration without noise influence.

和X^i(t)计算第i次迭代中无噪声影响的第t个快拍的接收 信号的近似值q^i(t)：(3) Utilization

The approximation of the received signal q ^i(t) of the t-th snapshot without noise influence in the ith iteration is calculated by using Xi ^(t) as follows:

Here, _gs ^(i-1)(t) represents the input function of the t snapshot of the i-1th iteration.

(4) Substitute γ ⁱ , (σ ² ) ⁱ and qi ^(t) into formula (28) to update the input function g _s ^{i (t)} of the t snapshot of the i-th iteration.

和输出阻尼系数

此外，再引入一个中间参数s，s的表达式为：(5) In order to ensure the convergence of GAMP iteration, the input damping coefficient is introduced

and output damping factor

In addition, an intermediate parameter s is introduced, and the expression of s is:

Among them, si ^(t) represents the value of the i-th iteration, which has no actual physical meaning and is just a mathematical equation. s ^{(i -1)(t)} represents the value of the i-1-th iteration.

(6) Use the updated input function g _si ^(t) to calculate the approximate value of the extended incident signal ri ^(t) and its noise variance of the t-th snapshot

Among them, A(θ)_sparse ^T represents the transpose of A(θ)_sparse, S ^T represents the transpose of S, and (g _s ^i(t) )' represents the derivative of the input function g _s ^i(t) .

(7) Substitute γ ⁱ , ri ^(t) and τ _ri ^(t) into formula (29) to update the output function g _x ^{i (t)} of the t snapshot of the i-th iteration.

(8) Update the variance and mean of the expanded incident signal according to the updated output function _gxi ^(t) :

表示第i次迭代输出的扩展后的入射信号的方差，(g_x ^i(t))'表 示输出函数g_x ^i(t)的导数，X^(i+1)(t)表示第i次迭代输出的扩展后的入射信号的 均值。in,

represents the variance of the expanded incident signal output at the i-th iteration, (g _x ^i(t) )' represents the derivative of the output function g _x ^i(t) , and X ^(i+1)(t) represents the mean of the expanded incident signal output at the i-th iteration.

进行GAMP迭代收敛 判断，GAMP迭代收敛条件为：(9) Using the GAMP iterative convergence condition, X ^(i+1)(t) and

Perform GAMP iteration convergence judgment. The GAMP iteration convergence condition is:

或迭代次数i达到GAMP算法的最 大迭代次数G，其中，ε_gamp表示GAMP算法迭代的归一化公差参数。

Or the number of iterations i reaches the maximum number of iterations G of the GAMP algorithm, where ε _gamp represents the normalized tolerance parameter of the GAMP algorithm iteration.

和X^(i+1)(t)，不满足 GAMP迭代收敛条件时，进入下一步。When the GAMP iteration convergence condition is met, the iteration ends and the output is

When X ^{(i+1)(t) and X(i+1)(t)} do not meet the GAMP iteration convergence condition, proceed to the next step.

分别计算第i+1次迭代中的先验误差γⁱ⁺¹和噪声方差 (σ²)ⁱ⁺¹，计算公式如下：(10) Update the iteration parameters based on the EM algorithm (expectation maximization algorithm), using the output of the i-th iteration X ^(i+1)(t) and

The prior error γ ⁱ⁺¹ and the noise variance (σ ² ) ⁱ⁺¹ in the i+1th iteration are calculated respectively, and the calculation formulas are as follows:

(σ ² ) ⁱ⁺¹ =||y ^(t) -A(θ)_sparse·X ^(i+1)(t) || ² (37)

进行EM迭代收敛判断， EM迭代收敛条件为：(11) Using the EM iterative convergence condition, X ^(i+1)(t) and

Perform EM iteration convergence judgment. The EM iteration convergence condition is:

或迭代次数i达到EM算法的最大迭 代次数K，其中，ε_em表示EM算法迭代的归一化公差参数。

Or the number of iterations i reaches the maximum number of iterations K of the EM algorithm, where ε _em represents the normalized tolerance parameter of the EM algorithm iteration.

不满足EM迭代收敛条件时，继续迭代。When the EM iteration convergence condition is met, the i+1th iteration is performed, and only steps (1) to (8) are performed in the i+1th iteration. The iteration is terminated and the mean value X ^(i+2)(t) and variance of the expanded incident signal of the i+1th iteration are output.

When the EM iteration convergence condition is not met, continue the iteration.

In step 5, each row of the restored expanded incident signal contains only W non-zero values, and the positions of the non-zero values correspond to the DOA directions of the incident signal x.

The present invention provides several sets of simulation experiments to verify the effect of the method of the present invention:

Simulation experiment 1:

A uniform linear array is used to receive the incident signal. The number of array elements of the receiving array is M=12. There are 4 signal sources, and the angles of the waves are -30.15, -10.23, 40.74, and 28.39. In the MMV model, the number of sampling snapshots T is 100. The method of the present invention and the existing TMSBL algorithm are used to estimate the direction of arrival of the signal. The normalized power spectra of the method of the present invention and the existing TMSBL algorithm in different noise environments are shown in Figures 3 and 4. It can be seen from the figure that when the signal-to-noise ratio SNR=10dB is low, the TMSBL algorithm that also takes time-related factors into consideration cannot estimate, while the method of the present invention can still estimate the direction of arrival of the signal more accurately; when the signal-to-noise ratio is high, that is, when SNR=40dB, the peak of the power spectrum of the TMSBL algorithm is still inaccurate, and there are false peaks. It can be seen that the method of the present invention is suitable for DOA model estimation in non-ideal environments.

Simulation experiment 2:

Assume that the number of array elements of the receiving array is M=10, there are two signal sources, the directions of the incoming waves are -11.21 and 0.74 respectively, and the number of snapshots T is 40. The predetermined grid point range is [-90°, 90°], and the resolution is 1°. The performance comparison of RMSE with the change of signal-to-noise ratio (SNR) is studied by using the method of the present invention and the existing TMSBL algorithm. For each set of simulation parameter settings, the root mean square error is obtained by averaging 500 Monte Carlo experiments.

FIG5 is a comparison diagram of the method of the present invention and the TMSBL algorithm on RMES. It can be seen that, taking into account the time correlation factor, the error of the method of the present invention is much smaller than that of the TMSBL algorithm at a low signal-to-noise ratio. Therefore, the accuracy of the method of the present invention is higher in non-ideal environments.

Simulation experiment 3:

While ensuring that other simulation conditions are the same as those of simulation experiment 2, the signal-to-noise ratio SNR is set to 30 dB, and the performance comparison of the CPU time with the change of the number of snapshots is studied by using the method of the present invention and the TMSBL algorithm. The results are shown in FIG6 . Compared with TMSBL, in the case of large snapshots, the time required by the method of the present invention is much less than that of TMSBL. The computational complexity of the method of the present invention is lower when the number of snapshots is large, and the computational efficiency is higher.

In summary, the method of the present invention solves the time correlation problem of the signal source in the MMV direction of arrival estimation under the Bayesian framework, uses GAMP to replace the matrix inversion part in the traditional SBL algorithm, reduces the computational complexity, and can estimate DOA more quickly and accurately. Considering non-ideal environments such as time correlation, it has the advantages of high resolution, low complexity, high accuracy, and good robustness.

The above description is only a preferred embodiment of the present invention. It should be pointed out that, for ordinary technicians in this technical field, several improvements and modifications can be made without departing from the technical principles of the present invention. These improvements and modifications should also be regarded as the protection scope of the present invention.

Claims (8)

1. A low-complexity signal direction of arrival estimation method based on sparse Bayes, characterized in that it comprises the following steps: Expanding the incident signal based on spatial domain grid partitioning; The Markov probability prior model is used to obtain the prior of the expanded incident signal; The input function and output function of the GAMP algorithm are obtained according to the prior of the expanded incident signal; Using the input function and output function of the GAMP algorithm, GAMP iteration is performed on the snapshot of the extended incident signal to obtain the restored extended incident signal; Obtaining the signal arrival direction according to the restored and expanded incident signal; The specific steps of incident signal expansion are as follows: Assume that the incident signal is x, and the receiving array consists of M uniform linear array antennas. When the receiving array receives T snapshots of the incident signal, the array output model is obtained: y＝A(θ)x+e Where y is the received signal of the receiving array, A(θ) is the array flow matrix, and e is the Gaussian white noise received by the receiving array; The spatial domain of the receiving array is gridded using the exhaustive method to obtain an overcomplete angle set.

in,

represents the nth supercomplete angle, H is the supercomplete angle set

The number of angles in, n = 1, 2, ..., H; based on

and A(θ) to obtain the expanded array flow matrix:

Among them, A(θ)_sparse represents the expanded array flow matrix; The overcomplete array output model is obtained based on the expanded array flow matrix:

in,

represents the received signal after expansion,

represents the incident signal after expansion; The steps to obtain the prior of the expanded incident signal are as follows: The expanded incident signal is modeled using the Markov probability prior model to obtain the relationship between the nth expanded incident signal at the tth snapshot and the t-1th snapshot:

in,

represents the nth expanded incident signal of the tth snapshot,

represents the nth expanded incident signal of the t-1th snapshot,

express

and

, β is the time correlation coefficient, β∈(-1,1), γ _n is the prior error of the incident signal after the nth expansion, t＝1,2,…,T; According to the relationship between the previous and next snapshots of the expanded incident signal, the influence of the previous snapshot is obtained.

And the impact of the next snapshot

in,

Represents the incident signal

The prior probability of

represents the mean of the forward pass of snapshot t-1 at the time of snapshot t of the nth extended incident signal,

represents the variance of the forward transmission of the t-1 snapshot at the t snapshot of the n-th expanded incident signal,

Represents the incident signal

The prior probability of

represents the mean value of the backward transmission of the t+1 snapshot of the incident signal after the n-th expansion,

represents the variance of the incident signal t+1 snapshot transmitted backwards at the time of snapshot t after the nth expansion; according to

and

Get the expanded incident signal

Priors:

2. According to a low-complexity signal direction-of-arrival estimation method based on sparse Bayesian according to claim 1, it is characterized in that the spatial domain of the receiving array is [-90°, 90°]. 3. The low-complexity signal direction of arrival estimation method based on sparse Bayes according to claim 1, characterized in that the nth expanded incident signal of the tth snapshot

The expression is as follows:

in,

4. The low-complexity signal direction of arrival estimation method based on sparse Bayes according to claim 1, characterized in that the input function of the GAMP algorithm is expressed as follows:

Where _gs represents the input function, q ^(t) represents the approximate value of the received signal of the t-th snapshot without noise influence,

represents the noise variance of q ^(t) , y ^(t) represents the received signal of the t-th snapshot receiving matrix, and ^σ2 represents the noise variance of y ^(t) ; The expression of the output function of the GAMP algorithm is as follows:

Among them, g _x represents the output function,

represents the approximate value of the n-th expanded incident signal of the t-th snapshot,

express

The noise variance. 5. The low-complexity signal direction of arrival estimation method based on sparse Bayes according to claim 4 is characterized in that the specific operation of obtaining the restored extended incident signal is as follows: Initialize the expanded incident signal of the t-th snapshot

Using the input function and output function of the GAMP algorithm

Perform GAMP iteration to update the mean and variance of the expanded incident signal; In each iteration process, the EM algorithm is used to update the iteration parameters, and the iteration convergence condition is used to judge the iteration convergence; Obtaining a restored extended incident signal according to the mean and variance of the extended incident signal that satisfy the iterative convergence judgment; Wherein, the iterative convergence condition includes the GAMP iterative convergence condition and the EM iterative convergence condition. 6. According to the low-complexity signal direction of arrival estimation method based on sparse Bayesian algorithm in claim 5, it is characterized in that, assuming that the maximum number of iterations of the GAMP algorithm is G, the iterative process of each iteration is specifically as follows: According to the output of the i-1th iteration, the extended incident signal of the tth snapshot in the i-th iteration is obtained.

Variance

mean Xi ^(t) , prior error γ ⁱ and noise variance of received signal (σ ² ) ⁱ , where i∈[1,G]; Let S = |A(θ)_sparse| ² , and use

The approximation qi ^(t) of the received signal of the t-th snapshot without noise influence in the ith iteration is calculated by using Xi ^(t) :

Where g _s ^(i-1)(t) represents the input function of the t snapshot of the i-1th iteration; Update the input function g _s ^i(t) of the t snapshot of the i-th iteration using γ ⁱ , (σ ² ) ⁱ and qi ^(t) ; Use the updated input function g _si ^(t) to calculate the approximate value ri ^(t) of the extended incident signal and its noise variance of the t-th snapshot.

Among them, A(θ)_sparse ^T represents the transpose of A(θ)_sparse,

is the input damping coefficient,

s ^(i-1)(t) represents the value of the i-1th iteration, S ^T represents the transpose of S, and (g _s ^i(t) )' represents the derivative of the input function g _s ^i(t) ; Using γ ⁱ , ri ^(t) and

Update the output function _gxi ^(t) of the t snapshot of the i-th iteration; Update the variance and mean of the expanded incident signal according to the updated output function _gxi ^(t) :

in,

represents the variance of the expanded incident signal output at the i-th iteration, (g _x ^i(t) )' represents the derivative of the output function g _x ^i(t) , X ^(i+1)(t) represents the mean of the expanded incident signal output at the i-th iteration,

is the output damping coefficient,

Using the GAMP iterative convergence condition, we can calculate X ^(i+1)(t) and

Perform GAMP iteration convergence judgment. When the GAMP iteration convergence condition is met, the iteration is terminated and the mean and variance of the incident signal after expansion of the i-th iteration are output. When the GAMP iteration convergence condition is not met, proceed to the next step. Based on the EM algorithm, using the output of the i-th iteration X ^(i+1)(t) and

The prior error γ ⁱ⁺¹ and the noise variance (σ ² ) ⁱ⁺¹ in the i+1th iteration are calculated respectively, and the calculation formulas are as follows:

(σ ² ) ⁱ⁺¹ =||y ^(t) -A(θ)_sparse·X ^(i+1)(t) || ² Using the EM iteration convergence condition, we can calculate X ^(i+1)(t) and

Perform EM iteration convergence judgment. When the EM iteration convergence condition is met, output the mean and variance of the expanded incident signal of the i+1th iteration, and end the iteration. When the EM iteration convergence condition is not met, continue the iteration. 7. The low-complexity signal direction of arrival estimation method based on sparse Bayes according to claim 6, characterized in that the GAMP iteration convergence condition is:

Or the number of iterations i reaches the maximum number of iterations G of the GAMP algorithm, where ε _gamp represents the normalized tolerance parameter of the GAMP algorithm iteration. 8. The low-complexity signal direction of arrival estimation method based on sparse Bayesian according to claim 6, wherein the EM iteration convergence condition is:

Or the number of iterations i reaches the maximum number of iterations K of the EM algorithm, where ε _em represents the normalized tolerance parameter of the EM algorithm iteration.

Images