CN110082708B - 非均匀阵列设计和波达方向估计方法 - Google Patents

非均匀阵列设计和波达方向估计方法 Download PDF

Info

Publication number
CN110082708B
CN110082708B CN201910139417.9A CN201910139417A CN110082708B CN 110082708 B CN110082708 B CN 110082708B CN 201910139417 A CN201910139417 A CN 201910139417A CN 110082708 B CN110082708 B CN 110082708B
Authority
CN
China
Prior art keywords
array
matrix
uniform
virtual
elements
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.)
Active
Application number
CN201910139417.9A
Other languages
English (en)
Other versions
CN110082708A (zh
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.)
Xidian University
Shaanxi University of Technology
Original Assignee
Xidian University
Shaanxi University of Technology
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 Xidian University, Shaanxi University of Technology filed Critical Xidian University
Priority to CN201910139417.9A priority Critical patent/CN110082708B/zh
Publication of CN110082708A publication Critical patent/CN110082708A/zh
Application granted granted Critical
Publication of CN110082708B publication Critical patent/CN110082708B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G01MEASURING; TESTING
    • G01SRADIO 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/00Direction-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/02Direction-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/14Systems for determining direction or deviation from predetermined direction
    • YGENERAL 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
    • Y02TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02DCLIMATE 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/00Reducing energy consumption in communication networks
    • Y02D30/70Reducing energy consumption in communication networks in wireless communication networks

Abstract

本发明公开了一种非均匀阵列设计和波达方向估计方法,主要解决现有技术布阵不灵活,计算复杂度较高的问题;其实现过程是:根据嵌套阵列最大自由度确定非均匀阵列首、尾位置系数;计算满足非均匀阵列差合矩阵
Figure DSA0000179553110000011
中包含所有虚拟阵元位置的位置系数;根据
Figure DSA0000179553110000012
最终得到非均匀阵列阵元位置;根据接收数据计算数据协方差矩阵并矢量化得到虚拟差分合成阵列接收数据r;对r去冗余、排序得到虚拟阵列接收数据
Figure DSA0000179553110000013
进而得到其满秩矩阵v;构造线性算子,通过估计线性算子得到信号子空间;构造选择矩阵得到旋转矩阵,最后通过旋转矩阵估计波达方向;本发明在同等条件下具有阵列配置灵活,无需对数据协方差矩阵进行特征分解,也不需对整个空域角度进行谱峰搜索的优点。

Description

非均匀阵列设计和波达方向估计方法
技术领域
本发明属于阵列信号处理技术领域,尤其涉及一种非均匀阵列设计和波达方向估计方法,通过灵活布阵,对多个空域目标源进行波达方向估计。
背景技术
在现代阵列信号处理领域中,波达方向估计(以下简称DOA估计)占有一个重要地位,传统DOA估计主要通过MUSIC方法、ESPRIT方法及其衍生方法进行。对于一个总阵元数为N的均匀线性阵列,采用传统DOA估计方法所能识别的最大空域目标个数为N-1个。在实际中,空域目标个数大于阵元数情况时有发生,采用少于目标个数的均匀线性阵列的传统DOA估计方法将失效。因此,空域目标个数大于阵元数的DOA估计被许多学者广泛研究,最近提出的嵌套阵列和互质阵列通过合理布阵构造一等效的虚拟阵列,提高阵列自由度(Degreeof freedom,DOF),通过利用虚拟阵列而不是原始阵列来进行DOA估计。嵌套阵列是由两个或多个具有不同阵元间隔的线性子阵组成,由它构造的虚拟阵列是一个完全填充的均匀线性阵列,但其第一个子阵阵元间距通常较小(阵元排布密度较大),会引起阵元间的互耦问题。互质阵列是由两个阵元数分别为M、N的均匀线性子阵组成,其中M、N为互质整数,虽然互质阵列能获得多余阵元数的自由度DOF,但其构造的虚拟阵列不是一个完全填充的均匀线性阵列,导致在DOA估计时,只能利用互质阵列构造的虚拟阵列的一部分,这就带来了自由度和阵列孔径的降低。
由于嵌套阵列、互质阵列每个子阵均由均匀线性子阵组成,这就导致在机载雷达或物理空间受限的情况下,时常出现无法为其子阵阵元安装找到合理位置的问题,同时关于非均匀阵列提出的新算法,发明名称为基于嵌套式子阵阵列的波达方向估计方法,公开(公告)号:CN105824002A,提出的对数据协方差矩阵特征分解后利用MUSIC算法谱峰搜索进行DOA估计;以及发明名称为嵌套阵列基于K-R子空间的波达方向估计方法,公开(公告)号:CN107544051A,提出的稀疏信号重构法利用了信号谱的稀疏性来进行DOA估计;这两种方法在进行DOA估计时,通常需要巨大的计算量。本发明的目的在于克服上述已有非均匀阵列以及DOA估计方法不足,设计一种新的非均匀阵列和波达方向估计方法,该非均匀阵列无需由均匀子阵列组成,且其生成的虚拟阵列是一个完全填充的均匀线性阵列,同时为降低计算量,通过构造一线性算子得到信号子空间进行DOA估计。本发明在保证可识别的空域目标个数大于阵元数前提下,提高了布阵的灵活性,降低DOA估计运算量。
发明内容
为实现上述目的,本发明的技术思路是:根据嵌套阵列经过最优布阵获取阵列最大自由度(最多虚拟阵元个数)思想,来确定非均匀阵列第一个阵元和最后一个阵元位置系数,通过计算非均匀阵列其余阵元位置系数,构造非均匀阵列;对非均匀阵列构造的虚拟差分合成阵列去冗余、排序得到虚拟阵列;构造线性算子,根据虚拟阵列进行DOA估计。其实现方案包括如下:
1)构造非均匀阵列:
1a)根据嵌套阵列最大自由度DOF确定非均匀阵列位置系数x1和xN
假定非均匀阵列阵元位置D=[d1,d2,...,dn,...,dN]=d×[x1,x2,...,xn,...,xN],根据嵌套阵列构造的虚拟阵列DULA=[-(N2(N1+1)-1)d,...,0,...,(N2(N1+1)-1)d],获取最大自由度DOF=2N2(N1+1)-1,得到非均匀阵列位置系数
Figure BSA0000179553130000021
则本发明阵列阵元位置为
Figure BSA0000179553130000022
其中d为阵元间隔,取值为入射信号最小半波长;N1与N2为嵌套阵列子阵阵元个数,N1=N2=N/2;xn表示非均匀阵列第n个阵元的位置系数,n=1,2,...,N,N为总阵元数;
Figure BSA00001795531300000210
元素为{1~(xN-1)}的随机递增互异整数;
1b)计算满足非均匀阵列差合矩阵
Figure BSA0000179553130000025
中包含所有虚拟阵元位置的位置系数
Figure BSA0000179553130000026
为保证本发明阵列在差合处理之后包含所有的虚拟阵元位置,且与嵌套阵列自由度一致,构造一向量
Figure BSA0000179553130000027
P=[D,D,...,D],向量P由N个D按行排列组成,再构造另一矩阵
Figure BSA0000179553130000028
即矩阵
Figure BSA0000179553130000029
为非均匀阵列产生的差合矩阵:
Figure BSA0000179553130000031
由于矩阵
Figure BSA0000179553130000032
为反对称矩阵,故只需研究其上三角形元素即可,将其按行排布得到一个行向量
Figure BSA0000179553130000033
为使非均匀阵列位置
Figure BSA0000179553130000034
所产生的虚拟差合阵列是完备的,即包含所有虚拟阵元位置,构造的差合矩阵
Figure BSA0000179553130000035
所对应的
Figure BSA0000179553130000036
中应包含
Figure BSA0000179553130000037
个不同元素;
1c)根据位置系数
Figure BSA0000179553130000038
获取非均匀阵列阵元位置;
在1a)中,随机从{1~(xN-1)}中选取N-2个互异递增数构成
Figure BSA0000179553130000039
得到一组均匀阵列阵元位置D,进行迭代,当D中阵元位置满足1b)中的差合矩阵
Figure BSA00001795531300000310
所对应
Figure BSA00001795531300000311
中包含
Figure BSA00001795531300000312
个不同值即可停止迭代,否则继续从1a)选取一组阵元位置,进行迭代,直至满足迭代停止条件,退出迭代,最终得到阵元位置
Figure BSA00001795531300000313
2)根据所设计非均匀阵列得到接收数据X,进而得到数据协方差矩阵Rx,向量化该协方差矩阵,得到虚拟差分合成阵列接收数据r;
3)根据虚拟差分合成阵列接收数据r,对r进行去冗余、排序操作,最终得到虚拟阵列接收数据
Figure BSA00001795531300000314
4)根据虚拟阵列接收数据
Figure BSA00001795531300000315
构造选择矩阵Jz
Figure BSA00001795531300000316
进行秩恢复操作得到满秩矩阵V;
5)根据虚拟阵列的满秩矩阵V,构造一线性算子Q,通过估计线性算子Q得到信号子空间ES
6)基于旋转因子不变法的思想定义两个选择矩阵Jg1和Jg2,获得ES1和ES2,通过ES1和ES2的相位关系获得旋转矩阵Φ,最后通过旋转矩阵估计波达方向
Figure BSA00001795531300000317
本发明与现有阵列结构相比较有如下优点:
1.本发明阵列虽与嵌套阵列自由度相同,但比嵌套布阵、互质布阵方式更灵活,不需要由两个或多个均匀线性子阵构成,且可产生完全填充型均匀线性虚拟阵列进行DOA估计,更适于工程实现;
2.本发明可满足在机载或者物理空间受限的环境下,从本发明中选取一种性能相同但却适合现场位置的非均匀阵列,进行灵活稀疏布阵,得到多余阵元数的自由度,扩大了阵列孔径,使得可估计的空域目标数远超过实际阵元数,具有更好的测角性能;
3.本发明无需对数据协方差矩阵进行特征分解,也不需通过对整个空域角度进行谱峰搜索,通过构造线性算子Q得到信号子空间ES,利用旋转因子不变法思想便可得到信号的DOA估计值,解决了在采样点数增多、空域网格步长划分较细,使得特征分解以及普峰搜索运算量大、耗时长、实时性差的问题;
4.本发明阵列为非均匀阵列,与传统均匀阵列相比,在阵列面积不变的情况下,能够有效减少阵元数目,降低成本,且阵元完全相同,更易于生产和降低成本。
附图说明
图1为本发明的三种六阵元非均匀阵列及虚拟阵列排布方式图;
图2为六阵元二阶嵌套阵列及产生虚拟阵元阵列排布方式图;
图3为本发明方法实现流程图;
图4为本发明方法空域目标DOA估计散点图;
图5为本发明和二阶嵌套阵、均匀线阵对DOA估计均方根误差与信噪比关系比较图;
具体实施方式
为了让本发明的上述和其它目的、特征及优点能更明显,下文特举本发明实施例,并配合所附图示,做详细说明如下:
参照图1为本发明三种六阵元非均匀阵列及产生虚拟阵列排布方式图,从图1可看到存在多种满足迭代条件的阵元位置系数
Figure BSA0000179553130000041
即存在多种阵元排布方式,要说明的是,嵌套阵列排布只是本发明阵列的一个子集;
参照图2为六阵元二阶嵌套阵列及产生虚拟阵元阵列排布方式图,本发明的非均匀阵列与嵌套阵列构造的虚拟阵列是一个完全填充的均匀线性阵列,即全部位置上虚拟阵元都是连续的,二者虚拟阵列的最大自由度DOF相同。因此,本发明阵列具有嵌套阵列的优点,不仅可识别多于阵元数的目标个数而且提高了目标方向识别的空间分辨力,同时具有布阵灵活的优点。
参照图3,由于嵌套阵列构造的虚拟阵列是一个完全填充的均匀线性阵列,可提高DOA估计的自由度DOF,本发明根据嵌套阵列经过最优布阵获取阵列最大自由度即最多虚拟阵元个数思想来设计一种新的非均匀阵列,构造一线性算子,利用连续的虚拟阵列进行DOA估计。为便于叙述,假定本发明阵列与二阶嵌套阵列总阵元数均为N,N为偶数,N为奇数的情况可类比得到,本发明方法步骤如下:
步骤一、根据嵌套阵列构造的虚拟阵列获取最大自由度DOF,计算非均匀阵列阵元位置系数,得到非均匀阵列;
1a)根据嵌套阵列最大自由度DOF确定非均匀阵列位置系数x1和xN
二阶嵌套阵为获取最大的阵列自由度,则每阶阵元数为N1=N2=N/2,其最大自由度为DOF=2N2(N1+1)-1,由其构造的虚拟阵列可等效为一个阵元位置为DULA=[-(N2(N1+1)-1)d,...,0,...,(N2(N1+1)-1)d]的均匀线阵,在总阵元数相同情况下为得到与二阶嵌套阵相同的阵列自由度,非均匀阵列位置矢量可表示为:D=[d1,d2,...,dn,...,dN]=d×[x1,x2,...,xn,...,xN];第一个阵元为参考阵元其位置系数x1=0,最后一个阵元位置系数
Figure BSA0000179553130000051
Figure BSA0000179553130000052
其中d为阵元间隔,取值为入射信号的最小半波长;xn表示非均匀阵列第n个阵元的位置系数,n=1,2,...,N,N为总阵元数;
Figure BSA0000179553130000053
Figure BSA00001795531300000513
元素为{1~(xN-1)}的随机递增互异整数;
1b)计算满足非均匀阵列差合矩阵
Figure BSA00001795531300000512
中包含所有虚拟阵元位置的位置系数
Figure BSA0000179553130000056
为保证本发明阵列与嵌套阵列构造的虚拟阵列一致,则本发明阵列必须满足在差合之后包含所有的虚拟阵元位置,即是完备的,其具体实现如下:
首先构造一向量
Figure BSA0000179553130000057
P=[D,D,...,D],构造另一个向量
Figure BSA0000179553130000058
Figure BSA0000179553130000059
即矩阵
Figure BSA00001795531300000510
为非均匀阵列产生的差合矩阵:
Figure BSA00001795531300000511
由于矩阵
Figure BSA0000179553130000061
为反对称矩阵,故只需研究其上三角形元素即可,将其按行排布得到一个行向量
Figure BSA0000179553130000062
为使非均匀阵列位置
Figure BSA0000179553130000063
所产生的虚拟差合阵列是完备的,即包含所有虚拟阵元位置,构造的差合矩阵
Figure BSA0000179553130000064
所对应的
Figure BSA0000179553130000065
中应包含
Figure BSA0000179553130000066
个不同元素;
1c)根据位置系数
Figure BSA0000179553130000067
获取非均匀阵列阵元位置
Figure BSA0000179553130000068
在1a)中,随机从{1~(xN-1)}中选取N-2个互异递增数构成
Figure BSA0000179553130000069
得到一组均匀阵列阵元位置D,进行迭代,当D中阵元位置满足1b)中的差合矩阵
Figure BSA00001795531300000610
所对应
Figure BSA00001795531300000611
中包含
Figure BSA00001795531300000612
个不同值即可停止迭代,否则继续从1a)选取一组阵元位置,进行迭代,直至满足迭代停止条件,退出迭代,最终得到阵元位置
Figure BSA00001795531300000613
步骤二、根据所设计的非均匀阵列,得到阵列接收数据X,进而得到数据协方差矩阵Rx,向量化该协方差矩阵,得到虚拟差分合成阵列接收数据r;
2a)本发明阵列接收数据信号为远场窄带信号,通过L次同步采样获取数据X;
2b)根据非均匀阵列接收数据X,估计得到数据协方差矩阵Rx
Figure BSA00001795531300000614
其中L表示快拍数,X∈CN×L,Rx∈CN×N,(·)H表示共轭转置;
2c)根据阵列协方差矩阵Rx,计算虚拟差分合成阵列接收数据r:
r=vec(Rx)
其中,
Figure BSA00001795531300000615
vec(·)表示对矩阵进行向量化操作;
步骤三、根据虚拟差分合成阵列接收数据r,对r进行去冗余、排序操作,最终得到虚拟阵列接收数据
Figure BSA00001795531300000616
由于接收数据r中包含噪声,故无法直接对r进行去冗余、排序操作,只有通过合理设计才能得到虚拟阵列接收数据
Figure BSA00001795531300000617
其原理是由于接收数据r元素位置与阵列产生的差合阵元位置一一对应,通过对差合阵元位置进行处理得到索引集,由索引集对r内部元素进行选取便可达到对r去冗余、排序目的,进而得到
Figure BSA00001795531300000618
具体实现如下:
3a)根据1)得到的非均匀阵列,构造一列向量
Figure BSA00001795531300000619
列向量
Figure BSA00001795531300000620
由1)中差合矩阵
Figure BSA00001795531300000621
的元素按行叠加组成一个列向量:
Figure BSA0000179553130000071
再构造另一向量U,其内部
Figure BSA0000179553130000072
个元素均匀连续变化,与非均匀阵列产生的虚拟阵元排布相同:
Figure BSA0000179553130000073
对列向量
Figure BSA0000179553130000074
中每一个元素
Figure BSA0000179553130000075
与向量U的Uj进行比较,在
Figure BSA0000179553130000076
时,记录
Figure BSA0000179553130000077
Figure BSA0000179553130000078
中的索引值i,并令Uj=N2,要说明的是,在进行比较时,只要满足
Figure BSA0000179553130000079
便更新索引集Γ=Γ∪{j},最终得到索引集Γ,其内部包含
Figure BSA00001795531300000710
Figure BSA00001795531300000711
个虚拟阵元的索引值,所对应的元素构成一个新的向量
Figure BSA00001795531300000712
其中,i=1,2,...,N2
Figure BSA00001795531300000713
N为总阵元数;
3b)根据对
Figure BSA00001795531300000714
排序得到新的索引集
Figure BSA00001795531300000715
获取虚拟阵列接收数据
Figure BSA00001795531300000716
由于索引集Γ所对应
Figure BSA00001795531300000717
中元素并不是按从小到大顺序排列的,因此通过对
Figure BSA00001795531300000718
排序即
Figure BSA00001795531300000719
向量
Figure BSA00001795531300000720
Figure BSA00001795531300000721
所对应F中的
Figure BSA00001795531300000722
个元素一一对应;
至此,得到排序后的索引值
Figure BSA00001795531300000723
通过接收数据r元素位置与阵列产生的差合阵元位置一一对应关系,用索引值
Figure BSA00001795531300000724
选取接收数据r便可达到去冗余、排序目的,最终得到虚拟阵列接收数据
Figure BSA00001795531300000725
其中,sort(·)表示排序运算,F表示内部元素按顺序排列后的向量,
Figure BSA00001795531300000726
表示F中每个元素对应索引值的集合;
步骤四、根据虚拟阵列接收数据
Figure BSA00001795531300000727
构造选择矩阵Jz
Figure BSA00001795531300000728
进行秩恢复操作得到满秩矩阵V,其具体过程如下:
4a)根据虚拟阵列接收数据
Figure BSA00001795531300000729
构造选择矩阵Jz
Figure BSA00001795531300000730
分割为
Figure BSA00001795531300000731
个子矩阵:
Figure BSA00001795531300000732
每个子矩阵用
Figure BSA00001795531300000733
表示,
Figure BSA00001795531300000734
由接收数据
Figure BSA00001795531300000735
Figure BSA00001795531300000736
Figure BSA00001795531300000737
行组成,其中Jz表示为:
Figure BSA00001795531300000738
4b)根据
Figure BSA00001795531300000739
叠加得到满秩矩阵V,V可表示为:
Figure BSA00001795531300000740
其中,
Figure BSA00001795531300000741
步骤五、非均匀阵列接收数据快拍数L、维度增加时,为了降低由特征分解获得信号子空间的运算量,根据满秩矩阵V,构造一线性算子,通过估计线性算子Q得到信号子空间ES,其具体过程如下:
5a)列出线性算子Q与信号子空间ES关系式:
Figure BSA0000179553130000081
其中,IK表示K×K的单位矩阵,K为目标个数,(·)H表示共轭转置,因此只需估计出线性算子即可获得信号子空间ES
5b)根据满秩矩阵V,估计线性算子Q,通过线性算子Q与信号子空间ES关系式获得信号子空间ES
为得到线性算子
Figure BSA0000179553130000082
分割满秩矩阵V:
Figure BSA0000179553130000083
其中
Figure BSA0000179553130000084
由V的前K行组成,
Figure BSA0000179553130000085
由V的
Figure BSA0000179553130000086
行组成,则线性算子能被估计如下:
Figure BSA0000179553130000087
获得线性算子Q后,通过
Figure BSA0000179553130000088
便可得到信号子空间;
步骤六、基于旋转因子不变法的思想定义两个选择矩阵Jg1和Jg2,获得ES1和ES2,通过ES1和ES2的相位关系获得旋转矩阵Φ,最后利用旋转矩阵估计波达方向
Figure BSA0000179553130000089
其具体实现如下:
6a)定义两个选择矩阵:
Figure BSA00001795531300000810
Figure BSA00001795531300000811
其中
Figure BSA00001795531300000812
表示
Figure BSA00001795531300000813
的单位矩阵,
Figure BSA00001795531300000814
表示
Figure BSA00001795531300000815
的全1向量,
Figure BSA00001795531300000816
Figure BSA00001795531300000817
分别表示
Figure BSA00001795531300000818
Figure BSA00001795531300000819
零矩阵;
6b)通过选择矩阵Jg1和Jg2获取信号子空间ES1和ES2,ES1=Jg1ES,ES2=Jg2ES,基于旋转因子不变法的思想,利用信号子空间ES1和ES2关系获得
Figure BSA00001795531300000820
对Ψ∈CK×K进行特征分解[v,q]=eig(Ψ)得到旋转矩阵Φ∈CK×K对角线元素:
Figure BSA00001795531300000821
其中eig(·)表示特征分解,
Figure BSA0000179553130000091
表示对矩阵求伪逆,diag(·)表示对角阵,v和q=[q1,q2,...,qk,...,qK]分别表示矩阵Ψ的特征向量矩阵和特征值矩阵,θk表示第k个信号的波达方向,λ为入射信号波长;
6c)通过旋转矩阵Φ估计波达方向
Figure BSA0000179553130000092
Figure BSA0000179553130000093
其中
Figure BSA0000179553130000094
表示k个波达方向估计值,Φk,k表示矩阵Φ的第k行第k列的元素,d为入射信号半波长,arg(·)表示取相位运算。
本发明的效果可以通过以下的仿真结果进一步说明:
仿真1:利用本发明非均匀阵列能估计数量大于接收阵列阵元数的目标数,图4仿真为6阵元非均匀阵列,阵元位置矢量D=d×[0,1,3,4,9,11],阵元个数为N=6,d/λ=1/2,对空间K=11个远场窄带目标进行方向识别,这11个目标的波大方向为(-50°,-40°,-30°,-20°,-10°,0°,10°,20°,30°40°,50°),在信噪比SNR=-10dB,快拍数为L=1024,进行50次Monte-Carlo实验;
从图4中可以看到,在低信噪比下,本发明非均匀阵列在阵元数为6的情况下能正确无误的识别出11个空域目标源,而传统六阵元均匀线阵在此情况下将不能进行目标方向识别;
仿真2:本发明非均匀阵列与二阶嵌套阵、均匀线阵在不同信噪比条件下的性能比较如图5所示,本发明非均匀阵列阵元个数为N=6,阵元位置矢量D=d×[0,1,3,4,9,11];二阶嵌套阵总阵元数为N1+N2=6,N1和N2分别为每阶子阵阵元个数,N1=N2=3,其二阶嵌套阵阵元位置矢量为Z=d×[0,1,2,3,7,11];均匀线阵总阵元数为6,,对空间K=2个远场窄带目标进行方向识别,两个目标的波达方向为(θ1,θ2)=(20°,40°),将目标DOA估计的结果和实际值比较,最后得到误差统计。实验信噪比SNR从0dB到40dB,快拍数为L=1024,进行500次Monte-Carlo实验;
从图5中可以看到,当信噪比在0dB下,本发明的非均匀阵列和二阶嵌套阵DOA估计精度在0.05°附近,均匀线阵的估计精度在0.34°附近,可以得到本发明方法优于均匀线阵;当信噪比在[0,20]dB时,本发明方法的估计精度优于二阶嵌套阵方法,随着信噪比的升高,均匀线阵估计性能上升,曲线下降加快,在40dB时三者性能趋近相同;
以上所述,仅是本发明的较佳实施例而已,并非对本发明做任何形式上的限制,虽然本发明已以较佳实施例揭露如上,然而并非用以限定本发明,任何熟悉本专业的技术人员,在不脱离本发明技术方案范围内,当可利用上述揭示的技术内容做出些许更动或修饰为等同变化的等效实施例,但凡是未脱离本发明技术方案的内容,依据本发明的技术实质对以上实施例所作的任何简单修改、等同变化与修饰,均仍属于本发明技术方案的范围内。

Claims (6)

1.一种非均匀阵列设计和波达方向估计方法,包含以下步骤:
1)构造非均匀阵列:
1a)根据嵌套阵列最大自由度DOF,确定非均匀阵列位置系数x1和xN
假定非均匀阵列阵元位置为:
D=[d1,d2,...,dn,...,dN]=d×[x1,x2,...,xn,...,xN]
根据嵌套阵列构造的虚拟阵列
DULA=[-(N2(N1+1)-1)d,...,0,...,(N2(N1+1)-1)d]
获取最大自由度DOF=2N2(N1+1)-1,得到非均匀阵列位置系数
Figure FSB0000204024750000011
则阵列阵元位置为
Figure FSB0000204024750000012
其中d为阵元间隔,取值为入射信号的最小半波长;N1与N2为嵌套阵列子阵阵元个数,N1=N2=N/2;xn表示非均匀阵列第n个阵元的位置系数,n=1,2,...,N,N为总阵元数;
Figure FSB0000204024750000013
由{1~(xN-1)}的随机递增互异整数组成;
1b)计算非均匀阵列差合矩阵
Figure FSB0000204024750000014
中包含所有虚拟阵元位置的位置系数
Figure FSB0000204024750000015
为保证阵列在差合处理之后包含所有的虚拟阵元位置,且与嵌套阵列自由度一致,构造一向量
Figure FSB0000204024750000016
P=[D,D,...,D],向量P由N个D按行排列组成,再构造另一矩阵
Figure FSB0000204024750000017
即矩阵
Figure FSB0000204024750000018
为非均匀阵列产生的差合矩阵:
Figure FSB0000204024750000019
Figure FSB0000204024750000021
由于矩阵
Figure FSB0000204024750000022
为反对称矩阵,故只研究其上三角元素即可,将其按行排布得到一个行向量
Figure FSB0000204024750000023
如下:
Figure FSB0000204024750000024
为使非均匀阵列位置
Figure FSB0000204024750000025
所产生的虚拟差合阵列是完备的,即包含所有虚拟阵元位置,构造的差合矩阵
Figure FSB0000204024750000026
所对应的
Figure FSB0000204024750000027
中应包含
Figure FSB0000204024750000028
个不同元素;
1c)根据位置系数
Figure FSB0000204024750000029
获取非均匀阵列阵元位置;
在1a)中,随机从{1~(xN-1)}中选取N-2个互异递增数构成
Figure FSB00002040247500000210
得到一组均匀阵列阵元位置D,进行迭代,当D中阵元位置满足1b)中的差合矩阵
Figure FSB00002040247500000211
所对应
Figure FSB00002040247500000212
中包含
Figure FSB00002040247500000213
个不同值即可停止迭代,否则继续从1a)选取一组阵元位置,进行迭代,直至满足迭代停止条件,退出迭代,最终得到阵元位置
Figure FSB00002040247500000214
2)根据所设计非均匀阵列得到接收数据X,进而得到数据协方差矩阵Rx,向量化该协方差矩阵,得到虚拟差分合成阵列接收数据r;
3)根据虚拟差分合成阵列接收数据r,对r进行去冗余、排序操作,最终得到虚拟阵列接收数据
Figure FSB00002040247500000215
4)根据虚拟阵列接收数据
Figure FSB0000204024750000031
构造选择矩阵Jz
Figure FSB0000204024750000032
进行秩恢复操作得到满秩矩阵V;
5)根据虚拟阵列的满秩矩阵V,构造一线性算子Q,通过估计线性算子Q得到信号子空间ES
6)基于旋转因子不变法的原理定义两个选择矩阵Jg1和Jg2,获得ES1和ES2,通过ES1和ES2的相位关系获得旋转矩阵Φ,最后通过旋转矩阵估计波达方向
Figure FSB0000204024750000033
2.根据权利要求1所述的非均匀阵列设计和波达方向估计方法,其中步骤2)中的虚拟差分合成阵列接收数据r,按如下步骤计算:
2.1)根据非均匀阵列接收数据X,估计得到数据协方差矩阵Rx
Figure FSB0000204024750000034
其中,L表示快拍数,X∈CN×L,Rx∈CN×N,(·)H表示共轭转置;
2.2)根据数据协方差矩阵Rx,计算虚拟差分合成阵列接收数据r∈N2×1:
r=vec(Rx)
其中vec(·)表示对矩阵进行向量化操作。
3.根据权利要求1所述的非均匀阵列设计和波达方向估计方法,其中步骤3)中的虚拟阵列数据接收矩阵
Figure FSB0000204024750000035
按如下步骤计算:
接收数据r中包含噪声,无法直接对r进行去冗余、排序操作,通过对差合阵元位置进行处理得到索引集,由索引集对r内部元素进行选取便可达到对r去冗余、排序目的,进而得到
Figure FSB0000204024750000036
具体实现如下:
3.1)根据1)得到的非均匀阵列,构造一列向量
Figure FSB0000204024750000041
列向量
Figure FSB0000204024750000042
由1)中差合矩阵
Figure FSB0000204024750000043
的元素按行叠加组成列向量
Figure FSB0000204024750000044
Figure FSB0000204024750000045
再构造另一向量U,其内部
Figure FSB0000204024750000046
个元素均匀连续变化,与非均匀阵列产生的虚拟阵元排布相同:
Figure FSB0000204024750000047
对列向量
Figure FSB0000204024750000048
中每一个元素
Figure FSB0000204024750000049
与向量U的Uj进行比较,在
Figure FSB00002040247500000410
时,记录
Figure FSB00002040247500000411
Figure FSB00002040247500000412
中的索引值i,并令Uj=N2,在进行比较时,只要满足
Figure FSB00002040247500000413
便更新索引集Γ=Γ∪{i},最终得到索引集Γ,其内部包含
Figure FSB00002040247500000414
Figure FSB00002040247500000415
个虚拟阵元的索引值,所对应的元素构成一个新的向量
Figure FSB00002040247500000416
其中,i=1,2,...,N2
Figure FSB00002040247500000417
Figure FSB00002040247500000418
N为总阵元数;
3.2)根据对
Figure FSB00002040247500000419
排序得到新的索引集
Figure FSB00002040247500000420
获取虚拟阵列接收数据
Figure FSB00002040247500000421
由于索引集Γ所对应
Figure FSB00002040247500000422
中元素并不是按从小到大顺序排列的,因此通过对
Figure FSB00002040247500000423
排序即
Figure FSB00002040247500000424
向量
Figure FSB00002040247500000425
Figure FSB00002040247500000426
所对应F中的
Figure FSB00002040247500000427
个元素一一对应;
至此,得到排序后的索引值
Figure FSB00002040247500000428
通过接收数据r元素位置与阵列产生的差合阵元位置一一对应关系,用索引值
Figure FSB00002040247500000429
选取接收数据r便可达到去冗余、排序目的,最终得到虚拟阵列接收数据
Figure FSB00002040247500000430
其中,sort(·)表示排序运算,F表示内部元素按顺序排列后的向量,
Figure FSB0000204024750000051
表示F中每个元素对应索引值的集合。
4.根据权利要求1所述的非均匀阵列设计和波达方向估计方法,其中步骤4)得到满秩矩阵V,按如下步骤计算:
4.1)根据虚拟阵列接收数据
Figure FSB0000204024750000052
构造选择矩阵Jz
Figure FSB0000204024750000053
分割为
Figure FSB0000204024750000054
个子矩阵:
Figure FSB0000204024750000055
每个子矩阵用
Figure FSB0000204024750000056
表示,
Figure FSB0000204024750000057
由接收数据
Figure FSB0000204024750000058
Figure FSB0000204024750000059
Figure FSB00002040247500000510
行组成,其中Jz表示为:
Figure FSB00002040247500000511
4.2)根据
Figure FSB00002040247500000512
叠加得到满秩矩阵V,V可表示为:
Figure FSB00002040247500000513
其中,
Figure FSB00002040247500000514
子阵序号
Figure FSB00002040247500000515
5.根据权利要求1所述的非均匀阵列设计和波达方向估计方法,其中步骤5)得到线性算子Q以及信号子空间ES,按如下步骤计算:
5.1)根据满秩矩阵V,列出线性算子Q与信号子空间ES关系式:
Figure FSB00002040247500000516
其中,IK表示K×K的单位矩阵,K为目标个数,(·)H表示共轭转置,因此只需估计出线性算子即可获得信号子空间ES
5.2)根据线性算子Q与信号子空间ES关系式,通过估计线性算子Q获得信号子空间ES
为得到线性算子
Figure FSB0000204024750000061
分割满秩矩阵V:
Figure FSB0000204024750000062
其中
Figure FSB0000204024750000063
由V的前K行组成,
Figure FSB0000204024750000064
由V的
Figure FSB0000204024750000065
行组成,则线性算子能被估计如下:
Figure FSB0000204024750000066
获得线性算子Q后,通过
Figure FSB0000204024750000067
便可得到信号子空间。
6.根据权利要求1所述的非均匀阵列设计和波达方向估计方法,其中步骤6)得到旋转矩阵Φ,最后利用旋转矩阵Φ估计波达方向
Figure FSB0000204024750000068
按如下步骤计算:
6.1)定义两个选择矩阵:
Figure FSB0000204024750000069
Figure FSB00002040247500000610
其中
Figure FSB00002040247500000611
表示
Figure FSB00002040247500000612
单位矩阵,
Figure FSB00002040247500000613
表示
Figure FSB00002040247500000614
的全1向量;
Figure FSB00002040247500000615
表示
Figure FSB00002040247500000616
Figure FSB00002040247500000617
零矩阵;
6.2)通过选择矩阵Jg1和Jg2获取信号子空间ES1和ES2,得到
Figure FSB0000204024750000071
对Ψ∈CK×K进行特征分解[v,q]=eig(Ψ)得到旋转矩阵Φ∈CK×K
Figure FSB0000204024750000072
其中,信号序号k=1,2,…,K,eig(·)表示特征分解,
Figure FSB0000204024750000073
表示对矩阵求伪逆,diag(·)表示矩阵对角化,v和q=[q1,q2,…,qk,…,qK]分别表示矩阵Ψ的特征向量矩阵和特征值矩阵,θk表示第k个信号的波达方向,λ为入射信号波长;
6.3)通过旋转矩阵Φ估计波达方向
Figure FSB0000204024750000074
Figure FSB0000204024750000075
其中
Figure FSB0000204024750000076
表示k个波达方向估计值,Φk,k表示矩阵Φ的第k行第k列元素,d为入射信号半波长,arg(·)表示取相位运算。
CN201910139417.9A 2019-02-25 2019-02-25 非均匀阵列设计和波达方向估计方法 Active CN110082708B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201910139417.9A CN110082708B (zh) 2019-02-25 2019-02-25 非均匀阵列设计和波达方向估计方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201910139417.9A CN110082708B (zh) 2019-02-25 2019-02-25 非均匀阵列设计和波达方向估计方法

Publications (2)

Publication Number Publication Date
CN110082708A CN110082708A (zh) 2019-08-02
CN110082708B true CN110082708B (zh) 2023-05-02

Family

ID=67413068

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201910139417.9A Active CN110082708B (zh) 2019-02-25 2019-02-25 非均匀阵列设计和波达方向估计方法

Country Status (1)

Country Link
CN (1) CN110082708B (zh)

Families Citing this family (11)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN110726967B (zh) * 2019-10-25 2021-08-03 北京理工大学 用于一维阵列测向的单边稀疏嵌套阵设计方法
CN110736959B (zh) * 2019-10-25 2021-07-09 北京理工大学 一种基于和差协同阵构建的平面互质阵列设计方法
CN111175691B (zh) * 2019-11-29 2021-11-05 北京理工大学 一种用于波达方向估计的双边稀疏嵌套阵设计方法
CN111751798A (zh) * 2020-07-22 2020-10-09 上海英恒电子有限公司 一种雷达测角方法
CN111896929B (zh) * 2020-08-28 2023-08-04 西安电子科技大学 非均匀mimo雷达的dod/doa估计算法
CN111983554A (zh) * 2020-08-28 2020-11-24 西安电子科技大学 非均匀l阵下的高精度二维doa估计
CN112327303A (zh) * 2020-10-22 2021-02-05 四川长虹电器股份有限公司 天线虚拟均匀线阵数据获取方法
CN112630724B (zh) * 2020-10-30 2022-11-08 哈尔滨工程大学 一种适用于uuv平台的高分辨目标方位估计方法
CN112632753B (zh) * 2020-12-01 2023-08-04 厦门大学 宽带非频变的非均匀间隔阵列天线方向图综合方法及装置
CN113435027A (zh) * 2021-06-23 2021-09-24 哈尔滨工程大学 一种高自由度低耦合稀疏线性阵列排布方法
CN113589224A (zh) * 2021-08-03 2021-11-02 宜宾电子科技大学研究院 一种基于增强嵌套阵的doa估计方法

Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105824002A (zh) * 2016-04-15 2016-08-03 西安电子科技大学 基于嵌套式子阵阵列的波达方向估计方法

Family Cites Families (12)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2002048853A (ja) * 2000-08-02 2002-02-15 Matsushita Electric Ind Co Ltd 電波到来方向推定装置及び指向性可変送受信装置
US20050195103A1 (en) * 2004-01-13 2005-09-08 Davis Dennis W. Phased arrays exploiting geometry phase and methods of creating such arrays
CN101272168B (zh) * 2007-03-23 2012-08-15 中国科学院声学研究所 一种信源数估计方法及其波达方向估计方法
CN101799535A (zh) * 2009-11-27 2010-08-11 西安电子科技大学 Mimo雷达目标方向的估计方法
CN102608565B (zh) * 2012-03-23 2014-01-08 哈尔滨工程大学 一种基于均匀圆阵列的波达方向估计方法
JP6004694B2 (ja) * 2012-03-26 2016-10-12 富士通テン株式会社 レーダ装置およびターゲット検出方法
JP6028388B2 (ja) * 2012-05-11 2016-11-16 富士通株式会社 到来方向推定装置、及び到来方向推定方法
CN103344939B (zh) * 2013-06-14 2015-10-28 西安交通大学 一种非相干及相干混合信号的二维波达方向估计方法
CN103399291B (zh) * 2013-07-22 2015-04-08 西安电子科技大学 基于快速稀疏恢复的超分辨波达方向估计方法
CN105403856B (zh) * 2015-10-30 2017-10-24 西安电子科技大学 基于嵌套式最小冗余阵列的波达方向估计方法
CN106019213B (zh) * 2016-05-09 2018-04-06 电子科技大学 一种部分稀疏l阵及其二维doa估计方法
CN108663653B (zh) * 2018-05-17 2020-04-07 西安电子科技大学 基于l形电磁矢量传感器阵列的波达方向估计方法

Patent Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105824002A (zh) * 2016-04-15 2016-08-03 西安电子科技大学 基于嵌套式子阵阵列的波达方向估计方法

Also Published As

Publication number Publication date
CN110082708A (zh) 2019-08-02

Similar Documents

Publication Publication Date Title
CN110082708B (zh) 非均匀阵列设计和波达方向估计方法
CN106324558B (zh) 基于互质阵列的宽带信号doa估计方法
CN110109050B (zh) 嵌套阵列下基于稀疏贝叶斯的未知互耦的doa估计方法
CN110244272B (zh) 基于秩一去噪模型的波达方向估计方法
CN111707985A (zh) 基于协方差矩阵重构的off-grid DOA估计方法
CN111983552B (zh) 一种基于差分共阵的嵌套阵列快速doa估计方法与装置
CN107037398B (zh) 一种二维music算法估计波达方向的并行计算方法
CN110837076A (zh) 一种基于张量分解的矢量水听器阵列方位估计方法
CN110895325B (zh) 基于增强四元数多重信号分类的到达角估计方法
CN111983554A (zh) 非均匀l阵下的高精度二维doa估计
CN111948599B (zh) 一种角度相关互耦影响下相干信号的高分辨率定位方法
CN113567913A (zh) 基于迭代重加权可降维的二维平面doa估计方法
CN112731280A (zh) 互质阵列混合噪声环境下的esprit-doa估计方法
CN113189538A (zh) 一种基于互质稀疏排列的三元阵列及其空间谱估计方法
CN112444773A (zh) 基于空域融合的压缩感知二维doa估计方法
CN111368256A (zh) 一种基于均匀圆阵的单快拍测向方法
CN113835063B (zh) 一种无人机阵列幅相误差与信号doa联合估计方法
CN111337872B (zh) 一种用于相干信源测向的广义doa矩阵方法
CN114648041A (zh) 一种基于平行稀疏阵列的二维欠定doa估计算法
CN111366891B (zh) 一种基于伪协方差矩阵的均匀圆阵单快拍测向方法
CN114609580A (zh) 一种基于非圆信号的无孔互质阵列设计方法
CN114325568A (zh) 脉冲噪声环境下基于bnc的嵌套阵列非圆信号doa估计方法
CN114460531A (zh) 一种均匀线阵music空间谱估计方法
CN113791379A (zh) 嵌套阵列非高斯环境下的正交匹配追踪doa估计方法
CN109061564B (zh) 基于高阶累积量的简化近场定位方法

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