一种基于零化约束的DOA与互耦的联合估计方法与流程

一种基于零化约束的doa与互耦的联合估计方法
技术领域
1.本发明属于信息技术领域，尤其是阵列信号处理技术领域，特别涉及存在阵元互耦情形下的阵列信号处理，具体涉及一种基于零化约束的doa与互耦的联合估计方法。


背景技术：

2.doa(波达方向定位技术)，是电子、通信、雷达、声呐等研究领域的行业内用语，通过处理接收到的回波信号，获取目标的距离信息和方位信息。doa估计属于阵列信号处理技术，是利用传感器阵列估计信号的波达方向，广泛应用于雷达、通信、声呐及医学诊断等军用和民用技术领域。基于参数建模的doa估计方法由于具有很高的分辨率而更受青睐。但这类方法的性能易受阵列的不完美性影响，例如阵元幅相误差，阵元互耦等。阵元互耦来源于波的传播特性，而难以避免。因此，研究如何校正互耦受到广泛关注。
3.为了校正互耦，公开文献1(z.ye and c.liu,“on the resiliency of music direction findingagainst antenna sensor coupling,”ieee transactions on antennas and propagation,vol.56,no.2,pp.371
–
380,2008)将阵列两端的阵元去除，以等效消除互耦的影响，只利用中间子阵列估计doa。该方法缺点在于没有充分利用阵列，存在孔径损失。公开文献2(f.sellone and a.serra,“a novel online mutual coupling compensation algorithm for uniform and linear arrays,”ieee transactions on signal processing,vol.55,no.2,pp.560
–
573,2007)采用交替最小化一个代价函数，以迭代计算的方法校正互耦。该方法缺点在于校正后的空间谱易存在伪峰，影响doa估计性能。此外，公开文献1和公开文献2的方法都需要互耦度的先验信息，当估计互耦度大于真实互耦度，两种方法的性能都会下降。


技术实现要素：

4.本发明的目的在于针对上述现有技术的不足，提出一种基于零化约束的doa与互耦的联合估计方法。
5.本发明方法包括以下步骤：
6.步骤(1)阵列布阵及互耦存在时的接收信号建模：
7.布设均匀线阵，阵元数量为m，相邻阵元间距为d；k个波长为λ，来自θ1,θ2,
…
,θk方向的窄带信号入射到均匀线阵，共接收到n个快拍，多块拍接收数据建模为：y＝cas e，y＝[y[1]，y[2]，
…
，y[n]]，信号矩阵s＝[s[1]，s[2]，
…
，s[n]]，噪声矩阵e＝[ε[1]，ε[2]，
…
，ε[n]]，c为互耦矩阵，c＝toeplitz([1,c1,c2,
…
,c
p
,0,
…
,0])，toeplitz(
·
)表示生成toeplitz矩阵，c1,c2,
…
,c
p
为互耦参数，p表示互耦度；
[0008]
其中，在第n个快拍时阵列的接收数据建模为：n＝1,2,
…
,n；其中，表示复数集合，流型矩阵a＝[a(θ1),a(θ2),
…
,a(θk)]，a(θ1),a(θ2),
…
,a(θk)表示导向矢量，第k个窄带信号导向矢量a(θk)的第m个元素为am(θk)＝exp[j(m-1)uk]，k＝1,2,
…
,k，m＝1,2,
…
,m，uk＝2πd cos(θk)/λ；s[n]为信号向量，ε[n]为噪声向
量。
[0009]
步骤(2)基于奇异值分解进行维度压缩：
[0010]
对y进行奇异值分解，则y经维度压缩后得到vs表示由y的k个最大奇异值对应的右奇异矢量构成的矩阵。
[0011]
步骤(3)建立基于零化约束的优化问题：
[0012]
优化问题表示为：其中，||
·
||2表示2范数，z的向量化z＝vec(z)，vec(
·
)表示将矩阵按列堆叠向量化，ik表示k阶单位矩阵，表示kronecker积，互耦矩阵c
′
＝toeplitz([1,c
t
,0,...,0])，互耦向量c＝[c1,c2,...,c
p
′
]
t
，p
′
表示互耦度的估计值，辅助参数ηk为γ的第m(k-1)至第mk个元素构成的向量，k＝1,2,...,k，待优化向量υ＝[c
t
,γ
t
,h
t
]
t
；多块拍矩阵；多块拍矩阵为构建toeplitz矩阵的算子，的第i行第j列元素为ηn的第i-j k 1个元素，i＝1,2,
…
,(m-k)，j＝1,2,
…
,(k 1)，零化滤波器ω为常数向量；s.t.表示约束条件，(
·
)h、(
·
)
t
分别表示取共轭转置和转置。
[0013]
步骤(4)采用迭代方法求解零化约束的优化问题：
[0014]
首先初始化其中c0＝0
p
′
，0
p
′
表示p
′×
1的零向量，γ的初始值γ0和h的初始值h0通过假设阵列无互耦并对z进行类tls-esptrit方法获得；然后进行迭代计算，在第q 1次迭代过程中得到υ
q 1
＝υq δυ，q≥0，其中为第q次迭代计算得到的结果，δυ为搜索向量，通过求解如下线性二次优化问题得到：
[0015]
其中，f向量第q次迭代得到的互耦矩阵g向量f向量jacobian矩阵g向量jacobian矩阵多块拍q矩阵即为γq的第m(k-1)至第mk个元素构成的向量，[
·
]
.,1:p
′
表示取矩阵的前p
′
列构成新的矩阵，toeplitz矩阵hankel矩阵
为的第1个元素，依此类推；为构建toeplitz矩阵的另一种算子，的第i行为i＝1,2,
…
,(m-k)，γ为反对角单位矩阵。
[0016]
步骤(5)doa与互耦估计：
[0017]
迭代收敛后，得到收敛解迭代收敛后，得到收敛解为互耦向量的估计值，为γ的估计值，为h的估计值，将带入得到多项式，求解该多项式的k个根得到x1,x2,...,xk，则doa估计为k＝1,2,...,k；其中，∠(
·
)取复数的幅角主值；根据得到估计的互耦矩阵
[0018]
优选地，步骤(3)中的ω取值为步骤(4)迭代计算中的h的初始值，即ω＝h0。
[0019]
本发明方法相对于现有技术具有以下有益效果：
[0020]
首先，本发明方法可利用全部的阵列数据，不存在孔径损失；其次，本发明方法无需互耦度的先验信息，即使互耦度被过估计，也不对性能有显著影响；第三，本发明方法的doa和互耦估计精度高。
附图说明
[0021]
图1为本发明方法的总流程图；
[0022]
图2为本发明方法与其他方法在不同互耦度估计值下doa估计性能比较示意图；
[0023]
图3为本发明方法与其他方法在不同信噪比下doa估计性能比较示意图；
[0024]
图4为本发明方法与其他方法在不同信噪比下互耦估计性能比较示意图；
[0025]
图5为本发明方法与其他方法在不同快拍数下doa估计性能比较示意图；
[0026]
图6为本发明方法与其他方法在不同快拍数下互耦估计性能比较示意图。
具体实施方式
[0027]
以下参照附图，对发明的具体技术方案和效果作进一步地详细说明。
[0028]
如图1所示，一种基于零化约束的doa与互耦的联合估计方法，具体步骤如下：
[0029]
步骤(1)阵列布阵及互耦存在时的接收信号建模：
[0030]
布设均匀线阵，阵元数量为m，相邻阵元间距为d；k个波长为λ，来自θ1,θ2,
…
,θk方向的窄带信号入射到均匀线阵，共接收到n个快拍，多块拍接收数据建模为：y＝cas e，y＝[y[1]，y[2]，
…
，y[n]]，信号矩阵s＝[s[1]，s[2]，
…
，s[n]]，噪声矩阵e＝[ε[1]，ε[2]，
…
，ε[n]]，c为互耦矩阵，对于均匀线阵，c具有对称toeplitz性质，且c＝toeplitz([1,c1,c2,
…
,c
p
,0,
…
,0])，toeplitz(
·
)表示生成toeplitz矩阵，c1,c2,
…
,c
p
为互耦参数，p表示互耦度。
[0031]
其中，在第n个快拍时阵列的接收数据建模为：n＝1,2,
…
,n；其中，表示复数集合，流型矩阵a＝[a(θ1),a(θ2),
…
,a(θk)]，a(θ1),a(θ2),
…
,a(θk)表示导向矢量，第k个窄带信号导向矢量a(θk)的第m个元素为am(θk)＝exp[j(m-1)uk]，k＝1,2,
…
,k，m＝1,2,
…
,m，uk＝2πd cos(θk)/λ；s[n]为信号向量，ε[n]为噪声向量，设定噪声为零均值高斯白噪声。
[0032]
步骤(2)基于奇异值分解进行维度压缩：
[0033]
对y进行奇异值分解，则y经维度压缩后得到vs表示由y的k个最大奇异值对应的右奇异矢量构成的矩阵。
[0034]
步骤(3)建立基于零化约束的优化问题：
[0035]
优化问题表示为：其中，||
·
||2表示2范数，z的向量化z＝vec(z)，vec(
·
)表示将矩阵按列堆叠向量化，ik表示k阶单位矩阵，表示kronecker积，互耦矩阵c
′
＝toeplitz([1,c
t
,0,...,0])，互耦向量c＝[c1,c2,...,c
p
′
]
t
，p
′
表示互耦度的估计值，辅助参数ηk为γ的第m(k-1)至第mk个元素构成的向量，k＝1,2,...,k，待优化向量υ＝[c
t
,γ
t
,h
t
]
t
；多块拍矩阵；多块拍矩阵为构建toeplitz矩阵的算子，的第i行第j列元素为ηn的第i-j k 1个元素，i＝1,2,
…
,(m-k)，j＝1,2,
…
,(k 1)，零化滤波器ω为常数向量，取值为迭代计算中h的初始值h0；s.t.表示约束条件，(
·
)h、(
·
)
t
分别表示取共轭转置和转置。
[0036]
步骤(4)采用迭代方法求解零化约束的优化问题；
[0037]
首先初始化其中c0＝0
p
′
，0
p
′
表示p
′×
1的零向量，γ的初始值γ0和h的初始值h0通过假设阵列无互耦并对z进行类tls-esptrit方法获得；然后进行迭代计算，在第q 1次迭代过程中得到υ
q 1
＝υq δυ，q≥0，其中为第q次迭代计算得到的结果，δυ为搜索向量，通过求解如下线性二次优化问题得到：
[0038]
其中，f向量第q次迭代得到的互耦矩阵g向量f向量jacobian矩阵g向量jacobian矩阵多块拍q矩阵即为γq的第m(k-1)至第mk个元素构成的向量，[
·
]
.,1:p
′
表示取矩阵的前p
′
列构成新的矩阵，分别为toeplitz矩阵和hankel矩阵，
为的第1个元素，依此类推；为构建toeplitz矩阵的另一种算子，的第i行为i＝1,2,
…
,(m-k)，γ为反对角单位矩阵。
[0039]
步骤(5)doa与互耦估计：
[0040]
迭代收敛后，得到收敛解迭代收敛后，得到收敛解为互耦向量的估计值，为γ的估计值，为h的估计值，将带入得到多项式，求解该多项式的k个根得到x1,x2,...,xk，则doa估计为k＝1,2,...,k；其中，∠(
·
)取复数的幅角主值；根据得到估计的互耦矩阵
[0041]
下面结合仿真实例对本发明方法的效果进行验证。
[0042]
仿真实例1：设置均匀线阵阵元数为m＝10，d＝λ/2，阵元互耦度p＝2，互耦参数为c1＝0.5 0.4j，c2＝0.1-0.03j。假设k＝2个入射信号，doa分别为θ1＝89
°
δθ1，θ2＝105
°
δθ2，δθ1与δθ2在[-0.5
°
,0.5
°
]内均匀分布。信噪比设置为15db，快拍数设置为100。将互耦度估计值从2扫描至9。比较本发明方法与背景技术中文献1和文献2的方法对doa估计的均方根误差上的性能差异。比较结果经500次蒙特卡洛试验平均后，在图2中显示。由于文献1方法的互耦度估计值最大只能取3，因此超过3的部分没有结果。由图2可看出，在互耦度过估计情况下，本发明方法的性能并未受显著影响，其doa估计的均方根误差仍小于θ1＝0.1
°
。其他两种方法则受较大影响。
[0043]
仿真实例2：设置快拍数为50个，信噪比从-10db扫描至15db，其余仿真条件与前述相同。比较各方法对doa和互耦的估计性能差异。为了提供性能参照，还提供了克拉美罗下限(crb)的结果。crb可通过文献：z.-m.liu and y.-y.zhou,“a unified framework and sparse bayesian perspective for direction-of-arrival estimation in the presence of array imperfections,”ieee transactions on signal processing,vol.61,no.15,pp.3786
–
3798,2013推得。仿真结果如图3和图4所示。从图中可以看出，对于doa和互耦的估计，本发明方法在各个信噪比的均方根误差达到最小，估计精度最高。并且在信噪比为5db时，本发明方法对doa的估计达到了crb。
[0044]
仿真实例3：设置信噪比为5db，快拍数从10扫描至100，其余仿真条件与前述相同。比较各方法对doa和互耦的估计性能差异。仿真结果如图5和图6所示。从图中可以看出，对于doa和互耦的估计，本发明方法在各个快拍数的均方根误差达到最小，估计精度最高。并且在快拍数为40时，本发明方法对doa的估计达到了crb。
[0045]
以上所述仅为本发名的较佳实施范例，并不用以限制本发明，凡在本发明的精神和原则之内所作的任何修改、等同替换和改进等，均应包含在本发明的保护范围之内。