ROC-MUSIC算法介绍

背景介绍

     在信号处理领域内,高斯分布假设一直占据着主导地位,因为其概率密度函数可以仅由均值和方差两个变量描述,这给信号的理论分析提供了极大的方便;而且基于高斯假设条件下的信号处理方法通常是线性的并且一般可以得到闭式最优解。然而,实际环境中存在许多非高斯特性的冲击噪声,如通信线路上的瞬间尖峰、由大气放电而引起的大气噪声以及雷达产生的部分杂波等。相比于高斯噪声,这类冲击噪声在时域上出现大量显著的尖峰脉冲特性,在概率密度上具有更加厚重的拖尾现象。Shao 和 Nikias 发现这种带有冲击性质的非高斯噪声可以用 α \alpha α稳定分布来描述。相对于高斯噪声,冲击噪声在实际中的用处更为广泛,因此关于冲击噪声背景下的研究也开始受到人们的关注,而DOA估计成为了研究的热点。
传统测向的方法较为突出的代表是MUSIC,ESPRIT算法等则在高斯噪声背景下表现出较好的性能,但是在冲击噪声背景下表现得往往不尽如人意,所以需要给出新的方法来实现DOA估计。

基于共变矩阵的ROC-MUSIC算法

     采用阵元数目为 M M M的均匀线阵对 N N N个信号进行测向,则阵元接收的信号可以表示为
[ x 1 ( t ) x 2 ( t ) ⋮ x M ( t ) ] = [ e − j ω 0 τ 11 e − j ω 0 τ 12 ⋯ e − j ω 0 τ 1 N e − j ω 0 τ 21 e − j ω 0 τ 22 ⋯ e − j ω 0 τ 2 N ⋮ ⋮ ⋮ e − j ω 0 τ M 1 e − j ω 0 τ M 2 ⋯ e − j ω 0 τ MN ] [ s 1 ( t ) s 2 ( t ) ⋮ s N ( t ) ] + [ n 1 ( t ) n 2 ( t ) ⋮ n M ( t ) ] \left[ \begin{matrix} {{x}_{1}}(t) \\ {{x}_{2}}(t) \\ \vdots \\ {{x}_{M}}(t) \\ \end{matrix} \right]=\left[ \begin{matrix} {{\text{e}}^{-\text{j}{{\omega }_{0}}{{\tau }_{11}}}} & {{\text{e}}^{-\text{j}{{\omega }_{0}}{{\tau }_{12}}}} & \cdots & {{\text{e}}^{-\text{j}{{\omega }_{0}}{{\tau }_{1N}}}} \\ {{\text{e}}^{-\text{j}{{\omega }_{0}}{{\tau }_{21}}}} & {{\text{e}}^{-\text{j}{{\omega }_{0}}{{\tau }_{22}}}} & \cdots & {{\text{e}}^{-\text{j}{{\omega }_{0}}{{\tau }_{2N}}}} \\ \vdots & \vdots & {} & \vdots \\ {{\text{e}}^{-\text{j}{{\omega }_{0}}{{\tau }_{M1}}}} & {{\text{e}}^{-\text{j}{{\omega }_{0}}{{\tau }_{M2}}}} & \cdots & {{\text{e}}^{-\text{j}{{\omega }_{0}}{{\tau }_{\text{MN}}}}} \\ \end{matrix} \right]\left[ \begin{matrix} {{s}_{1}}(t) \\ {{s}_{2}}(t) \\ \vdots \\ {{s}_{N}}(t) \\ \end{matrix} \right]+\left[ \begin{matrix} {{n}_{1}}(t) \\ {{n}_{2}}(t) \\ \vdots \\ {{n}_{M}}(t) \\ \end{matrix} \right] x1(t)x2(t)xM(t)=ejω0τ11ejω0τ21ejω0τM1ejω0τ12ejω0τ22ejω0τM2ejω0τ1Nejω0τ2Nejω0τMNs1(t)s2(t)sN(t)+n1(t)n2(t)nM(t)
将上式用矢量的形式表示则变为
X ( t ) = A ( θ ) S ( t ) + N ( t ) \mathbf{X}\left( t \right)=\mathbf{A}\left( \theta \right)\mathbf{S}\left( t \right)+\mathbf{N}\left( t \right) X(t)=A(θ)S(t)+N(t)
其中, X ( t ) \mathbf{X}\left( t \right) X(t)表示阵列接收到的 M × 1 M\times 1 M×1维的快拍数据, A ( θ ) \mathbf{A}\left( \theta \right) A(θ)表示 M × N M\times N M×N维的导向矢量, S ( t ) \mathbf{S}\left( t \right) S(t)表示 N × 1 N\times 1 N×1维的信号矢量, N ( t ) \mathbf{N}\left( t \right) N(t)表示 M × 1 M\times 1 M×1维的噪声矢量。
     假设接收到的信号和噪声都是服从 S α S S\alpha S SαS分布的随机变量,彼此之间独立,非相干。信号和噪声的共变矩阵分别为
{ Γ s = d i a g ( γ s 1 , γ s 2 , ⋯   , γ s N ) Γ N = d i a g ( γ n , γ n , ⋯   , γ n ) \left\{ \begin{aligned} & {{\mathbf{\Gamma }}_{s}}=diag\left( {{\gamma }_{{{s}_{1}}}},{{\gamma }_{{{s}_{2}}}},\cdots ,{{\gamma }_{{{s}_{N}}}} \right) \\ & {{\mathbf{\Gamma }}_{N}}=diag\left( {{\gamma }_{n}},{{\gamma }_{n}},\cdots ,{{\gamma }_{n}} \right) \\ \end{aligned} \right. {Γs=diag(γs1,γs2,,γsN)ΓN=diag(γn,γn,,γn)
     阵列的共变矩阵为
[ Γ i k ] = [ x i ( t ) , x k ( t ) ] α = ∑ j = 1 N a i ( θ j ) a k < α − 1 > ( θ j ) γ s j + γ n δ i k i , k = 1 , 2 ⋯   , M \left[ {{\mathbf{\Gamma }}_{ik}} \right]={{\left[ {{x}_{i}}(t),{{x}_{k}}(t) \right]}_{\alpha }}=\sum\limits_{j=1}^{N}{{{a}_{i}}}\left( {{\theta }_{j}} \right)a_{k}^{<\alpha -1>}\left( {{\theta }_{j}} \right){{\gamma }_{{{s}_{j}}}}+{{\gamma }_{n}}{{\delta }_{ik}}\quad i,k=1,2\cdots ,M [Γik]=[xi(t),xk(t)]α=j=1Nai(θj)ak<α1>(θj)γsj+γnδiki,k=1,2,M
     矢量表示形式为
Γ X = A Γ s A H + γ n I {{\mathbf{\Gamma }}_{X}}=\mathbf{A}{{\mathbf{\Gamma }}_{s}}{{\mathbf{A}}^{H}}+{{\gamma }_{n}}\mathbf{I} ΓX=AΓsAH+γnI
     根据共变矩阵的定义可以求出
Γ i k = E { x i ( t ) ∣ x k ( t ) ∣ p − 2 x k ∗ ( t ) } E { ∣ x k ( t ) ∣ p } 1 < p ≤ 2 {{\mathbf{\Gamma }}_{ik}}\text{=}\frac{E\left\{ {{x}_{i}}(t){{\left| {{x}_{k}}(t) \right|}^{p-2}}x_{k}^{*}(t) \right\}}{E\left\{ {{\left| {{x}_{k}}(t) \right|}^{p}} \right\}}\quad 1<p\le 2 Γik=E{xk(t)p}E{xi(t)xk(t)p2xk(t)}1<p2
通过上式可以看出,矩阵 Γ \mathbf{\Gamma } Γ M − N M-N MN个较小的特征值对应 γ n {{\gamma }_{n}} γn,其特征矢量张成的噪声子空间与矩阵 A \mathbf{A} A的列矢量正交,信号子空间则是由 N N N个较大特征值对应的特征矢量张成的空间组成。因此ROC-MUSIC算法的流程可以简述为以下过程
1)计算 M × M M\times M M×M维的矩阵 Γ \mathbf{\Gamma } Γ,其元素为
Γ i k = ∑ t = 1 K x i ( t ) ∣ x k ( t ) ∣ p − 2 x k ∗ ( t ) ∑ t = 1 K ∣ x k ( t ) ∣ p 1 < p < α ≤ 2 {{\mathbf{\Gamma }}_{ik}}\text{=}\frac{\sum\limits_{t=1}^{K}{{{x}_{i}}(t){{\left| {{x}_{k}}(t) \right|}^{p-2}}x_{k}^{*}(t)}}{\sum\limits_{t=1}^{K}{{{\left| {{x}_{k}}(t) \right|}^{p}}}}\quad 1<p<\alpha \le 2 Γik=t=1Kxk(t)pt=1Kxi(t)xk(t)p2xk(t)1<p<α2
2)对矩阵 Γ \mathbf{\Gamma } Γ进行奇异值分解,构成 M × ( M − N ) M\times \left( M-N \right) M×(MN)矩阵 E n = ( e ^ N + 1 , e ^ N + 2 , ⋯   , e ^ M ) {{E}_{n}}=\left( {{{\hat{e}}}_{N+1}},{{{\hat{e}}}_{N+2}},\cdots ,{{{\hat{e}}}_{M}} \right) En=(e^N+1,e^N+2,,e^M) ,其中 e ^ N + 1 , e ^ N + 2 , ⋯   , e ^ M {{\hat{e}}_{N+1}},{{\hat{e}}_{N+2}},\cdots ,{{\hat{e}}_{M}} e^N+1,e^N+2,,e^M Γ \mathbf{\Gamma } Γ中最小的 M − N M-N MN个奇异值对应的奇异向量
3)计算谱
S ( θ ) = 1 a H ( θ ) E n E n H a ( θ ) S(\theta )=\frac{1}{{{\mathbf{a}}^{\text{H}}}(\theta ){{\mathbf{E}}_{n}}\mathbf{E}_{n}^{\text{H}}\mathbf{a}(\theta )}\quad S(θ)=aH(θ)EnEnHa(θ)1
其中 a ( θ ) = [ 1 , e − j 2 π d 2 λ sin ⁡ θ , e − j 2 π d 3 λ sin ⁡ θ , ⋯   , e − j 2 π d M λ sin ⁡ θ ] T \mathbf{a}(\theta )={{\left[ 1,{{\text{e}}^{-\text{j}2\pi \frac{{{d}_{2}}}{\lambda }\sin \theta }},{{\text{e}}^{-\text{j}2\pi \frac{{{d}_{3}}}{\lambda }\sin \theta }},\cdots ,{{\text{e}}^{-\text{j}2\pi \frac{{{d}_{M}}}{\lambda }\sin \theta }} \right]}^{\text{T}}} a(θ)=[1,ej2πλd2sinθ,ej2πλd3sinθ,,ej2πλdMsinθ]T
4) S ( θ ) S\left( \theta \right) S(θ) N N N个局部极大值即为波达方向角 { θ 1 , θ 2 , ⋯   , θ N } \left\{ {{\theta }_{1}},{{\theta }_{2}},\cdots ,{{\theta }_{N}} \right\} {θ1,θ2,,θN}的估计。

仿真实验

仿真实验所用到的信号和噪声分别是远场窄带信号和冲击噪声,采用阵元数为10的均匀线阵,阵元间距为半波长,三个独立的入射信号,入射角度分别为-5°,15°和35°,广义信噪比为3dB,具体参数如下:

仿真参数参数值
阵元数目10
信源数目3
阵元间距半波长
快拍数1000
α \alpha α1.4
γ \gamma γ1

对于服从对称 α \alpha α稳定分布的冲击噪声来说,它的分布特性由特征指数 α \alpha α和尺度参数 γ \gamma γ 决定,但是当 α ∈ ( 0 , 2 ] \alpha \in (0,2] α(0,2]时,其方差是无界的,因而可以用信号的平均功率和冲击噪声的尺度参数的比值即广义信噪比(GSNR)代替一般的信噪比,其表达式如下
GSNR = 10 lg ⁡ ( 1 γ K ∑ t = 1 K ∣ s ( t ) 2 ∣ ) \text{GSNR}=10\lg \left( \frac{1}{\gamma K}\sum\limits_{t=1}^{K}{\left| s{{\left( t \right)}^{2}} \right|} \right) GSNR=10lg(γK1t=1Ks(t)2)
其中, K K K为快拍数。当 α = 2 \alpha \text{=}2 α=2时,上式退化为一般的信噪比。

冲击噪声的具体参考:冲击噪声介绍
实验结果如下:
在这里插入图片描述

     冲击噪声背景下,传统MUSIC算法不能对入射角度进行准确的估计,这是因为冲击噪声会对入射信号进行较大的干扰,从而对结果造成影响甚至是误判,而经过对噪声进行重新建模分析的ROC-MUSIC算法很好地区分出了入射信号的角度。
     令 α = 2 \alpha=2 α=2,此时的冲击噪声退化为高斯噪声,相同仿真条件下得到如下结果
在这里插入图片描述
     当冲击噪声退化成高斯噪声时,无论是传统的MUSIC算法还是ROC-MUSIC算法能较好地估计出入射信号的角度,但是此时传统MUSIC算法表现出较好性能,这也和前面的理论分析一致。

参考文献:

[1]P. Tsakalides and C. L. Nikias, “The robust covariation-based MUSIC (ROC-MUSIC) algorithm for bearing estimation in impulsive noise environments,” in IEEE Transactions on Signal Processing, vol. 44, no. 7, pp. 1623-1633, July 1996, doi: 10.1109/78.510611.
[2]赵大勇,刁鸣,杨丽丽,陈超.冲击噪声背景下的动态DOA跟踪[J].山东大学学报(工学版),2010,40(01):133-138.
[3]李洪升,张瑞丰,杜宇,王永孝.冲击噪声环境中波束域DOA估计方法研究[J].现代防御技术,2017,45(03):93-97.

  • 2
    点赞
  • 6
    收藏
    觉得还不错? 一键收藏
  • 9
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论 9
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值