CN110297209B - 一种基于平行互质阵列时空扩展的二维波达方向估计方法 - Google Patents
一种基于平行互质阵列时空扩展的二维波达方向估计方法 Download PDFInfo
- Publication number
- CN110297209B CN110297209B CN201910274526.1A CN201910274526A CN110297209B CN 110297209 B CN110297209 B CN 110297209B CN 201910274526 A CN201910274526 A CN 201910274526A CN 110297209 B CN110297209 B CN 110297209B
- Authority
- CN
- China
- Prior art keywords
- array
- data
- matrix
- signal
- 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
Links
- 239000011159 matrix material Substances 0.000 title claims abstract description 52
- 238000000034 method Methods 0.000 title claims abstract description 29
- 238000005070 sampling Methods 0.000 claims abstract description 11
- 238000003491 array Methods 0.000 claims abstract description 10
- 238000005314 correlation function Methods 0.000 claims abstract description 8
- 238000005457 optimization Methods 0.000 claims abstract description 5
- 239000013598 vector Substances 0.000 claims description 27
- 230000010363 phase shift Effects 0.000 claims description 3
- 230000001934 delay Effects 0.000 claims description 2
- 238000005311 autocorrelation function Methods 0.000 abstract description 6
- 230000003416 augmentation Effects 0.000 description 5
- 238000010586 diagram Methods 0.000 description 4
- 238000012545 processing Methods 0.000 description 4
- 238000000342 Monte Carlo simulation Methods 0.000 description 2
- 238000005516 engineering process Methods 0.000 description 2
- 238000011160 research Methods 0.000 description 2
- 241000584803 Xanthosia rotundifolia Species 0.000 description 1
- 230000003190 augmentative effect Effects 0.000 description 1
- 238000004891 communication Methods 0.000 description 1
- 230000021615 conjugation Effects 0.000 description 1
- 238000010276 construction Methods 0.000 description 1
- 230000007547 defect Effects 0.000 description 1
- 238000013461 design Methods 0.000 description 1
- 230000035800 maturation Effects 0.000 description 1
- 230000017105 transposition Effects 0.000 description 1
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/04—Details
- G01S3/10—Means for reducing or compensating for quadrantal, site, or like errors
-
- 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
- Physics & Mathematics (AREA)
- Engineering & Computer Science (AREA)
- General Physics & Mathematics (AREA)
- Radar, Positioning & Navigation (AREA)
- Remote Sensing (AREA)
- Measurement Of Velocity Or Position Using Acoustic Or Ultrasonic Waves (AREA)
Abstract
本发明提供一种基于平行互质阵列时空扩展的二维波达方向估计方法。该方法选取和增加特定的参考阵元,使互质阵列的子阵分别对参考阵元求互相关,得到不同时滞的相关函数,并将该自相关函数作为新的接收数据,然后由新的接收数据伪采样得到伪数据矩阵,对伪数据矩阵矢量化,将问题转化为稀疏优化问题求解得到信号的二维波达方向。本发明充分的利用了接收数据的时间‑空间信息,最大限度的拓展出虚拟阵列,具有更大的阵列孔径和更高的自由度,波达方向估计的性能更好。
Description
技术领域
本发明涉及阵列信号处理领域,尤其指一种基于平行互质阵列时空扩展的二维波达方向估计方法。
背景技术
阵列信号处理技术在众多领域已得到广泛应用,而阵列信号处理的基本问题之一是空间信号波达方向估计(DOA估计)。这一技术应用广泛,在通信、雷达等领域发挥着越来越重要的作用。在过去的几十年里人们不断提出了一系列高效的
DOA估计算法。近几年来,随着压缩感知理论的不断成熟和广泛应用,基于稀疏重构技术的DOA估计算法受到更多的关注和研究。使用均匀线阵估计信号的波达方向时,阵列的自由度直接受限于阵元个数,所能估计的信号源数目有限。由于非均匀阵列能够获得更大的阵列孔径和更高的自由度,并且阵元的摆放位置也具有更大的灵活性,所以研究非均匀阵列的设计和相应的DOA估计算法具有很好的实际意义和应用价值。互质并行阵列是一种非均匀阵列,由两个平行放置的互质子阵列组成,阵元间距为半波长的整数倍。通过矢量化子阵的互协方差矩阵,得到阵元数目更多的虚拟阵列,在实际物理阵列数一定的情况下,增大阵列的孔径,使得所能估计的信号源数目大于实际物理阵元数目,并且DOA估计的准确率也进一步提高。在文献An Efficient Dictionary Learning-Based 2-D DOAEstimationWithout Pair Matching for Co-Prime Parallel Arrays,IEEE Access,2018(以下简称文献[1]),提出了一种利用时间和空间信息的互质并行阵列扩展的方法,利用自相关函数的共轭对称性,构造出更多的虚拟阵元,进一步拓展了阵列孔径。
上述文献[1]中的方法虽然对互质阵列进行了虚拟扩展,但仅仅利用了部分空间信息来构造虚拟阵列,如果所有的时间-空间信息都能被利用,可以得到更多的虚拟阵元数目,那么阵列孔径和自由度将会进一步增加,DOA估计也会更加准确。
针对上述问题,我们提出了一个基于Co-prime并行阵列时空扩展的二维DOA估计方法。首先,我们在互质并行阵列的基础上增加了两个额外的阵元,利用信号自相关函数的共轭对称性,得到不同时滞的自相关函数,将该自相关函数作为新的接收数据,并求新的数据的协方差矩阵,最后通过稀疏重构的方法还原信号。所提出的方法可以利用空间和时间信息并形成共轭增强时空虚拟阵列,该方法能在物理阵元数目一定的情况下,利用空间和时间信息完全的扩展互质并行阵列,最大限度的提高阵列孔径和自由度以及得到更加准确的DOA估计值。
发明内容
本发明的目的在于克服现有方法存在的各种不足,首先提供一种基于平行互质阵列时空扩展的波达方向估计方法,在文献[1]原有互质阵列的基础上,增加两个物理阵元,充分利用接收信号的时间和空间信息对物理阵列进行虚拟扩展,形成共轭增强的时空虚拟阵列,最大限度的提高阵列孔径和自由度,得到更加准确的DOA估计值。
另外本发明的另一目的在于提供一种基于互质阵列扩展的阵列结构。
为了实现上述目的,本发明可以通过以下技术方案实现。
一种基于平行互质阵列时空扩展的二维波达方向估计方法,其包括:选取和增加参考阵元,使互质阵列的子阵分别对参考阵元求互相关,得到不同时滞的相关函数,并将该自相关函数作为新的接收数据,拓展出更多的虚拟阵元,由新的接收数据伪采样得到伪数据矩阵,然后对伪接收数据的协方差矩阵矢量化,最后通过稀疏重构的方法还原信号。
进一步地,基于互质阵列扩展的阵列结构包括一个平行互质阵列以及两个单独的阵元,互质并行阵列包含两个水平均匀线阵,均匀线阵分别有M1、M2个阵元,并且M1和M2互质;线阵1的阵元间距为线阵2的阵元间距为/>两个子阵的间距为/>以坐标系来说明各个阵元的相对位置,以线阵1第一个参考阵元为原点,线阵1位于坐标系水平轴上,线阵2位于线阵1上方(不仅限于上方,此仅举例说明本阵列的结构),间距为/>在互质并行阵列基础上增加的两个单独的阵元分别位于坐标点(λ,0),/>
进一步地,所述的一种基于平行互质阵列时空扩展的波达方向估计方法,其特征在于包含如下步骤:
步骤一:子阵1信号接收X1(t)=A1S(t)+Z1(t)
其中,A1=A1(α)=[a1(α1),a1(α2),…,a1(αK)]表示子阵1的流形矩阵,K表示信号源的个数,α=[α1,α2,…,αK],αk表示第k个信号与y轴之间的夹角,表示第k个信号的方向向量,S(t)=[s1(t),s2(t),…,sK(t)]T表示信号矢量,/>表示第k个信号源,Z1(t)表示子阵1接收到的零均值高斯白噪声矢量,方差为/>
步骤二:子阵2信号接收X2(t)=A2φS(t)+Z2(t),子阵2信号相对于子阵1信号有φ的相移,φ中包含俯仰角度信息;
其中,A2=A2(α)=[a2(α1),a2(α2),…,a2(αK)]表示子阵2的流形矩阵,K表示信号源的个数,α=[α1,α2,…,αK],αk表示第k个信号与y轴之间的夹角,表示第k个信号的方向向量,S(t)=[s1(t),s2(t),…,sK(t)]T表示信号矢量,/>表示第k个信号源,Z2(t)表示子阵2接收到的零均值高斯白噪声矢量,方差为/>
β=[β1,β2,…,βK],βk表示第k个信号与x轴之间的夹角。
步骤三:新增参考阵元信号接收,βi表示信号源俯仰角,
其中,y(t)表示阵元接收到的K个信号源数据之和,si(t)表示第i个信号源数据,i表示第i个信号源,z3(t)表示阵元/>接收到的噪声数据。w(t)表示阵元(λ,0)接收到的K个信号源数据之和,z4(t)表示阵元(λ,0)接收到的噪声数据。
由步骤一和步骤二得到的X1(t)和X2(t)分别对原有参考阵元(0,0),以及新增加的参考阵元/>(λ,0)的数据求互相关,以及利用相关函数的共轭对称性得到如下数据:
其中M1、M2为子阵1和子阵2的阵元数目,下角标M1表示阵列的第M1个阵元,下角标M1+1表示阵列的第M1+1个阵元,以此类推。Rxx(τ)表示子阵1和子阵2的数据分别对参考阵元(0,0),的数据作相关运算;Rxy(τ)表示子阵1的数据对参考阵元/>的数据作相关运算;Rxw(τ)表示子阵2的数据对参考阵元/>的数据作相关运算,字母τ表示相关函数的时延,求相关函数的公式定义为:
其中,表示第k个信号源的相关函数。(an,bn)为参考阵元的坐标,(am,bm)为子阵中的阵元坐标,m,n表示阵元的位置。
步骤四:构建虚拟阵列,由步骤三求得的数据,得到虚拟阵列接收信号:
对R1(τ)和R2(τ)以不同时延τ采样得到伪数据矩阵:
其中seqv=[Rs(Ts),Rs(2Ts),…,Rs(NpTs)],Np是虚拟矩阵接收数据的伪快拍数,Ts是虚拟矩阵接收数据的伪采样周期。
求r1和r2的互协方差矩阵:
RC=E[r2r1 H]
互协方差矩阵可以表示为:
其中,Reqv=E[seqv(seqv)H]是一个对角阵,该矩阵第k个元素为且φReqv也是对角阵。
矢量化协方差矩阵,则有:
其中,⊙表示矩阵的Khatri-Rao积,向量u包含对角矩阵φReqv中的对角元素。由于Reqv是实值对角矩阵,φ中对角元素包含角度β的数据,而向量u中的元素包含φ的信息,所以一旦求得向量u,就能计算出角度β。
稀疏重构算法求解。构造过完备字典Θ来替代A,将空间划分为多个网格{θ1,θ2,…,θD}(D>>K),真正的波达方向必定位于其中的网格附近。过完备字典由所有可能的方向{θ1,θ2,…,θD}(D>>K)构成。根据稀疏重构理论,将上述问题构建成一个稀疏优化问题,则有:
其中ρ=[ρ1,ρ2,…,ρD]T对应u中所有可能的方向。只要求解出ρ中的非零元素,就能求解出所有波达方向。具体来说,ρ中的非零元素所在的网格对应角度α,ρ中的非零元素的值对应于角度β,η是平衡稀疏度的正则化参数。
可以通过MATLAB软件的工具包求解得到ρ,ρ中非零元素的位置对应过完备基Θ中网格的位置,此即为角度α的估计值。ρ中非零元素值为u中元素的估计值,用如下公式求解β:
βk=arccos(-angle(uk)/π),k=1,2,…,K,
在本式中,uk为ρ中第k个非零元素的值,βk表示信号第k信号的俯仰角,K表示信号源的个数。
与现有技术相比,本发明具有如下优点:本发明与现有互质并行阵列相比,充分的利用时间-空间信息扩展了虚拟阵列,在物理阵元数相同的条件下,能拓展出比原有互质并行阵列更多的虚拟阵元,具有更大的阵列孔径与自由度,均方误差更小,具有更高的实用性。
附图说明
图1是现有技术中互质并行阵列结构图。
图2是实施例中增加两个阵元的互质并行阵列结构图。
图3是实施例中互质并行阵列时空拓展后的阵元结构图。
图4是实施例中增加两个阵元的互质并行阵列时空拓展结构图。
图5是实施例中不同信噪比下各个算法的DOA估计性能比较图。
图6是实施例中不同采样快拍数下各个算法的DOA估计性能比较图。
图7是本发明基于平行互质阵列时空扩展的二维波达方向估计方法的流程图。
具体实施方式
下面结合实施例子及附图对本发明的实施作进一步详细的描述,但本发明的实施方式不限于此,需指出的是,以下若有未特别详细说明之过程或符号,均是本领域技术人员可参照现有技术理解或实现的。
1、本实施例的阵列结构
本发明提供的基于互质阵列扩展的阵列结构如图2所示,该阵列结构包括一个平行互质阵列以及两个单独的阵元。互质并行阵列包含两个水平均匀线阵,均匀线阵分别有M1、M2个阵元,且M1与M2互质。线阵1的阵元间距为λ表示信号的波长,线阵2的阵元间距为/>两个子阵的间距为/>以坐标系来说明各个阵元的相对位置,线阵1位于坐标系水平轴上,线阵2位于线阵1上方(不仅限于上方,此仅举例说明本阵列的结构),间距为/>在互质并行阵列基础上增加的两个单独的阵元分别位于坐标点(λ,0),/>
2、基于平行互质阵列时空扩展的二维波达方向估计方法
本发明提供的基于互质并行阵列时空扩展的波达方向估计方法包含以下几个具体步骤:
步骤一:子阵1信号接收。K个窄带目标信号源分别为S1,S2,…,SK,信号源矢量记作s(t)=[s1(t)s2(t)…sK(t)]T。信号的波长为λ。这些目标声源对应于水平坐标轴的方向角分别为α=[α1,α2,…,αK],对应于垂直坐标轴的方向角分别为β=[β1,β2,…,βK]。αk表示第k个信号源的水平角,βk表示第k个信号源的俯仰角,则该阵元接收到的信号为:
其中,sk(t)表示第k个信号源,且有1≤k≤K。n1(t)表示第一个阵元上的噪声。
接收信号满足窄带条件,即当信号延迟远小于带宽倒数时,延迟作用相当于使基带信号产生一个相移。那么第m个阵元在同一时刻接收到的信号为:
其中λ表示信号源波长,αk表示第k个信号源的水平角,nm(t)表示第m个阵元上的噪声。将各阵元的接收信号排列成列向量形式,则整个子阵1的接收信号可用以下矢量式子表示:
X1(t)=A1S(t)+Z1(t) (1)
其中,为M1×K的导
向矢量矩阵,为M1×1的接收信号矩阵,S(t)=[s1(t),s2(t),…,sK(t)]T为K×1的源信号矩阵,sk(t)表示第k个信号源,且有1≤k≤K,Z1(t)为M1×1的噪声矩阵。
步骤二:子阵2信号接收,线阵2相对于线阵1有固定相位的延时,可以得到线阵2的信号接收模型为:
X2(t)=A2φS(t)+Z2(t) (2)
其中,为M2×K的导
向矢量矩阵,为M1×1的接收信号矩阵,S(t)=[s1(t),s2(t),…,sK(t)]T为K×1的源信号矩阵,sk(t)表示第k个信号源,且有1≤k≤K,Z2(t)为M1×1的噪声矩阵,/>
步骤三:共轭增广时空处理,用坐标点(am,bm),(an,bn)表示阵列中的两个不同的阵元,对应的接收信号表示为xm(t)和xn(t),τ表示时延,其协方差矩阵表示为:
其中,表示为:
表示为:
因此,公式(3)可以表示为:
(1)对阵元(0,0)共轭增广
使(an,bn)=(0,0),则公式(3)可以表示为:
令 对子阵1和子阵2分别有:
R(1)(τ)=A1Rs(τ),R(2)(τ)=A2φRs(τ) (5)
因为所以/>分别用/>和R(1-)(τ)表示/>和R(1)(-τ)*的最后M1-1行,则有:
(2)对阵元共轭增广
使则公式(3)可以表示为:
令
则有:R(1')(τ)=A1φ*Rs(τ),R(2')(τ)=A2Rs(τ) (8)
因为所以有/>
令R(1'-)(τ)表示(R(1')(-τ))*的后面M1-1行,和R(2'-)(τ)分别表示/>和(R(2')(-τ))*的后M2-1行。
(3)对阵元共轭增广,该阵元的接收信号为:
使则公式(3)可以表示为:
令则有:
R(3)(τ)=A1φRs(τ) (11)
(4)对阵元(λ,0)共轭增广,该阵元的接收信号为:
使(an,bn)=(λ,0),则公式(3)可以表示为:
令则有:
因为所以/>令:
步骤四,构建虚拟阵列,由公式(5)(6)(8)(9)(11)(12),则有:
R1(τ)和R2(τ)可以视为虚拟阵列的接收信号,由于虚拟阵列是通过不同时滞和空间位置的互相关获得的,因此虚拟阵列称为空间-时间虚拟阵列,如图4所示为时空拓展后的虚拟阵列结构。
如图3所示,在未添加两个单独阵元,仅对阵元(0,0)进行共轭增广的情况下,子阵1和子阵2分别扩展出两组虚拟阵元。如图4所示,在增加两个单独阵元,并对阵元(λ,0),/>进行共轭增广的情况下,子阵1和子阵2分别再次拓展出一组阵元。比较图4和图3,与文献[1]中互质并行阵列的时空拓展相比,本发明方法拓展的虚拟阵列包含更多传感器元素和更大的孔径,可以增大自由度。
由R1(τ)和R2(τ)构造伪数据矩阵:
其中seqv=[Rs(Ts),Rs(Ts),…,Rs(NpTs)],Np是伪快拍数,Ts是伪采样周期。
因此,求得r1和r2的互协方差矩阵:
互协方差矩阵可以表示为:
其中,Reqv=E[seqv(seqv)H],φReqv是对角阵。
矢量化协方差矩阵,则有:
其中,向量u包含对角矩阵φReqv中的对角元素。由于Reqv是实值矩阵,φ中对角元素包含角度β,向量u中的元素包含φ,所以一旦求得u,就能计算出角度β。
步骤五,稀疏重构算法求解,在本步骤中,我们使用稀疏重构算法来求解u。首先构造过完备字典Θ来替代A,将空间划分为多个网格{θ1,θ2,…,θD}(D>>K),真正的波达方向必定位于其中的网格中,那么过完备字典由所有可能的方向{θ1,θ2,…,θD}(D>>K)构成。根据稀疏重构理论,将上述问题看成稀疏优化问题,则有:
其中ρ=[ρ1,ρ2,…,ρD]T对应u中所有可能的方向。只要求解出ρ中的非零元素,就能求解出所有波达方向。具体来说,ρ中的非零元素所在的网格对应角度α,ρ中的非零元素的值对应于角度β。η是平衡稀疏度的正则化参数。
公式(19)可以通CVX工具包求解ρ,ρ中非零元素的位置对应过完备基Θ中网格的位置,此即为α的估计值。ρ中非零元素值为u的估计值,用如下公式求解β:
βk=arccos(-angle(uk)/π),k=1,2,…,K, (20)
uk为ρ中第K个非零元素的值。
以上算法流程图可以用图7表示。
本发明的主要工作步骤具体如下:
步骤1:信号源数目设定为K=5,信号源与水平轴夹角为25°,36°,47°,58°,69°,与竖直轴水平夹角为83°,71°,59°,47°,35°,子阵1的阵元数目M1为4,子阵2的阵元数目M2为3。
步骤2:子阵1接收到的K个信号源的数据X1(t)=A1S(t)+z1(t),子阵2接收到的K个信号源数据为X2(t)=A2φS(t)+N(t),两个单独阵元接收到的K个信号源数据为以下又分为四个子步骤:
(1)子阵1对阵元(0,0)求自相关得到R(1)(τ),对R(1)(τ)求共轭转置得到R(1-)(τ)。子阵2对阵元(0,0)求自相关得到R(2)(τ)。
(2)子阵1对阵元求自相关并共轭转置得到R(1'-)(τ),子阵2对阵元/>求自相关得到R(2')(τ),对R(2')(τ)求共轭转置得到R(2'-)(τ)。
(3)子阵1对求自相关得到R(3)(τ)。
(4)子阵2对(λ,0)求自相关得到R(4-)(τ)。
得到虚拟阵列的接收数据:
步骤3:对虚拟数据矩阵进行采样,得到r1=[R1(Ts),R1(Ts),…,R1(NpTs)],r2=[R2(Ts),R2(Ts),…,R2(NpTs)],其中伪快拍数设置为Np=500.Ts=1/4000s。求虚拟阵列接收数据的互协方差矢量化r=vec(Rc)。
步骤4:稀疏重构求解。构造过完备基Θ,设置网格为[0:0.2:180],则有
Θ=AA1*⊙AA1,求解稀疏优化问题可以用CVX工具包求解ρ,ρ中非零元素的位置对应过完备基Θ中网格的位置,此即为α的估计值。ρ中非零元素值为u的估计值,用如下公式求解β:
βk=arccos(-angle(uk)/π),k=1,2,…,K,
在本式中,uk为ρ中第k个非零元素的值,βk表示信号第k信号的俯仰角,K表示信号源的个数。
步骤5,性能比较。实施例中用均方根误差(RMSE)来评估各个算法的性能,其定义为
如图5所示,采样快拍数为500,参考信噪比SNR从-10dB到20dB变化,间隔为5dB,每个信噪比进行1000次蒙特卡洛实验。“Method in[1]M1=4M2=3”为文献[1]中的阵列结构下(如图1),阵元数M1=4,M2=3时的均方误差,“Method in[1]M1=4M2=5”为文献[1]中的阵列结构下(如图1),M1=4,M2=5时的均方误差,“Our Method M1=4M2=3”为本发明所提供的阵列结构及扩展方法下(如图2),增加两个单独阵元,且M1=4,M2=3时的均方误差。可以看出,本发明提供的算法性能优于文献[1]中算法,并且在阵元数目相等的情况下,性能同样优于原有算法。如图6所示,参考信噪比SNR为10dB,采样快拍数从100到1000变化,间隔为100,每个快拍数进行1000次蒙特卡洛实验,可以看出本发明所提供的算法同样性能较好。本发明在增加两个单独阵元的条件下,最大限度的对原有阵列进行了时间-空间拓展,得到更多的虚拟阵元,提高阵列孔径,具有更小的均方误差。
Claims (2)
1.一种基于平行互质阵列时空扩展的二维波达方向估计方法,其特征在于包括:选取和增加参考阵元,使互质阵列的子阵分别对参考阵元求互相关,得到不同时滞的相关函数,并将该相关函数作为新的接收数据,拓展出更多的虚拟阵元,由新的接收数据伪采样得到伪数据矩阵,然后对伪接收数据的协方差矩阵矢量化,最后通过稀疏重构的方法还原信号;
所述方法包含如下步骤:
步骤一:子阵1信号接收x1(t)=A1s(t)+z1(t)
其中,A1=A1(α)=[a1(α1),a1(α2),…,a1(αK)]表示子阵1的流形矩阵,K表示信号源的个数,α=[α1,α2,…,αK],αk表示第k个信号与y轴之间的夹角,表示第k个信号的方向向量,S(t)=[s1(t),s2(t),…,sK(t)]T表示信号矢量,/>表示第k个信号源,z1(t)表示子阵1接收到的零均值高斯白噪声矢量,方差为/>
步骤二:子阵2信号接收X2(t)=A2φS(t)+Z2(t),子阵2信号相对于子阵1信号有φ的相移,φ中包含俯仰角度信息;
其中,A2=A2(α)=[a2(α1),a2(α2),…,a2(αK)]表示子阵2的流形矩阵,K表示信号源的个数,α=[α1,α2,…,αK],αk表示第k个信号与y轴之间的夹角,表示第k个信号的方向向量,s(t)=[s1(t)s2(t)…sK(t)]T表示信号矢量,/>表示第k个信号源,z2(t)表示子阵2接收到的零均值高斯白噪声矢量,方差为σn 2;β=[β1,β2,…,βK],βk表示第k个信号与x轴之间的夹角;
步骤三:新增参考阵元信号接收
其中,y(t)表示阵元接收到的K个信号源数据之和,si(t)表示第i个信号源数据,i表示第i个信号源,z3(t)表示阵元/>接收到的噪声数据;w(t)表示阵元(λ,0)接收到的K个信号源数据之和,z4(t)表示阵元(λ,0)接收到的噪声数据;
由步骤一和步骤二得到的X1(t)和X2(t)分别对原有参考阵元(0,0),以及新增加的参考阵元/>(λ,0)的数据求互相关,以及利用相关函数的共轭对称性得到如下数据:
R(1-)(τ)=(R(1)(-τ))*
R(2'-)(τ)=(R(2')(-τ))*
R(3)(τ)=[Rx1y(τ),…,Rx1y(τ)]T
R(1'-)(τ)=(R(1')(-τ))*
R(4-)(τ)=(R(4-)(-τ))*
其中M1、M2为子阵1和子阵2的阵元数目,下角标M1表示阵列的第M1个阵元,下角标M1+1表示阵列的第M1+1个阵元,以此类推;Rxx(τ)求两个子阵阵元数据的相关函数,Rxy(τ)表示子阵1的数据对参考阵元的数据作相关运算;Rxw(τ)表示子阵2的数据对参考阵元的数据作相关运算,字母τ表示相关函数的时延,求相关函数的公式定义为:
其中,表示第k个信号源的相关函数;(an,bn)为参考阵元的坐标,(am,bm)为子阵中的阵元坐标,m,n表示阵元的位置;
步骤四:构建虚拟阵列,由步骤三求得的数据,得到虚拟阵列接收信号:
其中,
对R1(τ)和R2(τ)以不同时延τ采样得到伪数据矩阵:
其中seqv=[Rs(Ts),Rs(2Ts),…,Rs(NpTs)],Np是虚拟矩阵接收数据的伪快拍数,Ts是虚拟矩阵接收数据的伪采样周期,R1(Ts)表示对虚拟数据矩阵R1(τ)在Ts时刻的值,R2(Ts)表示对虚拟数据矩阵R2(τ)在Ts时刻的值,
求r1和r2的互协方差矩阵:
互协方差矩阵可以表示为:
其中,Reqv=E[seqv(seqv)H]是一个对角阵,且φReqv也是对角阵,β=[β1,β2,…,βK],βk表示第k个信号与x轴之间的夹角;
矢量化协方差矩阵,则有:
其中,⊙表示矩阵的Khatri-Rao积,向量u为对角矩阵φReqv中的对角元素;由于Reqv是实值对角矩阵,φ中对角元素包含角度β的数据,而向量u中的元素包含φ的信息,所以一旦求得向量u,就能计算出角度β;
稀疏重构算法求解:构造过完备字典Θ来替代A,将空间划分为多个网格{θ1,θ2,…,θD}(D>>K),真正的波达方向必定位于其中的网格附近;过完备字典由所有可能的方向{θ1,θ2,…,θD}(D>>K)构成;根据稀疏重构理论,将上述问题构建成一个稀疏优化问题,则有:
其中ρ=[ρ1,ρ2,…,ρD]T对应u中所有可能的方向;只要求解出ρ中的非零元素,就能求解出所有波达方向;具体来说,ρ中的非零元素所在的网格对应角度α,ρ中的非零元素的值对应于角度β,η是平衡稀疏度的正则化参数;
通过MATLAB软件的工具包求解得到ρ,ρ中非零元素的位置对应过完备基Θ中网格的位置,此即为角度α的估计值;ρ中非零元素值为u中元素的估计值,用如下公式求解β:
βk=arccos(-angle(uk)/π),k=1,2,…,K,
在本式中,uk为ρ中第k个非零元素的值,K表示信号源的个数。
2.根据权利要求1所述的一种基于平行互质阵列时空扩展的二维波达方向估计方法,其特征在于基于互质阵列扩展的阵列结构包括一个平行互质阵列以及两个单独的阵元,平行互质阵列包含两个水平均匀线阵,均匀线阵分别有M1、M2个阵元,并且M1和M2互质;线阵1的阵元间距为线阵2的阵元间距为/>两个子阵的间距为/>以坐标系来说明各个阵元的相对位置,以线阵1第一个参考阵元为原点,线阵1位于坐标系水平轴上,线阵2位于线阵1上方,间距为/>在平行互质阵列基础上增加的两个单独的阵元分别位于坐标点(λ,0),/>
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201910274526.1A CN110297209B (zh) | 2019-04-08 | 2019-04-08 | 一种基于平行互质阵列时空扩展的二维波达方向估计方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201910274526.1A CN110297209B (zh) | 2019-04-08 | 2019-04-08 | 一种基于平行互质阵列时空扩展的二维波达方向估计方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN110297209A CN110297209A (zh) | 2019-10-01 |
CN110297209B true CN110297209B (zh) | 2023-07-18 |
Family
ID=68026502
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201910274526.1A Active CN110297209B (zh) | 2019-04-08 | 2019-04-08 | 一种基于平行互质阵列时空扩展的二维波达方向估计方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN110297209B (zh) |
Families Citing this family (5)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN110736959B (zh) * | 2019-10-25 | 2021-07-09 | 北京理工大学 | 一种基于和差协同阵构建的平面互质阵列设计方法 |
CN112363107B (zh) * | 2020-09-21 | 2023-08-22 | 华南理工大学 | 基于互质阵列的混合信号波达方向估计方法 |
CN112763972B (zh) * | 2020-12-30 | 2023-05-09 | 长沙航空职业技术学院 | 基于稀疏表示的双平行线阵二维doa估计方法及计算设备 |
CN113325364A (zh) * | 2021-07-15 | 2021-08-31 | 金陵科技学院 | 一种基于数据压缩的空时联合测向方法 |
CN114325559B (zh) * | 2021-11-23 | 2023-03-28 | 电子科技大学 | 一种用于二维doa估计的互质平面阵列的布阵方法 |
Family Cites Families (9)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US7714782B2 (en) * | 2004-01-13 | 2010-05-11 | Dennis Willard Davis | Phase arrays exploiting geometry phase and methods of creating such arrays |
CN104749552A (zh) * | 2015-03-21 | 2015-07-01 | 西安电子科技大学 | 基于稀疏重构的互质阵列波达方向角估计方法 |
CN106324558B (zh) * | 2016-08-30 | 2021-04-02 | 东北大学秦皇岛分校 | 基于互质阵列的宽带信号doa估计方法 |
JP6853642B2 (ja) * | 2016-09-26 | 2021-03-31 | パナソニック株式会社 | レーダ装置 |
CN107015190A (zh) * | 2017-03-01 | 2017-08-04 | 浙江大学 | 基于虚拟阵列协方差矩阵稀疏重建的互质阵列波达方向估计方法 |
CN107315160B (zh) * | 2017-05-03 | 2019-06-18 | 浙江大学 | 基于内插虚拟阵列信号原子范数最小化的互质阵列波达方向估计方法 |
CN107329108B (zh) * | 2017-05-03 | 2019-07-23 | 浙江大学 | 基于内插虚拟阵列协方差矩阵Toeplitz化重建的互质阵列波达方向估计方法 |
CN108594164B (zh) * | 2017-11-30 | 2020-09-15 | 山东农业大学 | 一种平面阵列doa估计方法及设备 |
CN109581275B (zh) * | 2018-12-13 | 2023-09-29 | 华南理工大学 | 基于非圆信号和三维正交阵的二维水下doa估计方法和装置 |
-
2019
- 2019-04-08 CN CN201910274526.1A patent/CN110297209B/zh active Active
Also Published As
Publication number | Publication date |
---|---|
CN110297209A (zh) | 2019-10-01 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN110297209B (zh) | 一种基于平行互质阵列时空扩展的二维波达方向估计方法 | |
CN109932680B (zh) | 一种基于平移互质阵列的非圆信号波达方向估计方法 | |
Zheng et al. | DOA estimation for coprime linear arrays: An ambiguity-free method involving full DOFs | |
CN108872929B (zh) | 基于内插虚拟阵列协方差矩阵子空间旋转不变性的互质阵列波达方向估计方法 | |
Han et al. | An ESPRIT-like algorithm for coherent DOA estimation | |
CN109655799B (zh) | 基于iaa的协方差矩阵向量化的非均匀稀疏阵列测向方法 | |
CN107037392B (zh) | 一种基于压缩感知的自由度增加型互质阵列波达方向估计方法 | |
CN107561484B (zh) | 基于内插互质阵列协方差矩阵重建的波达方向估计方法 | |
CN106972882B (zh) | 基于虚拟域空间功率谱估计的互质阵列自适应波束成形方法 | |
CN109507636B (zh) | 基于虚拟域信号重构的波达方向估计方法 | |
CN107703478B (zh) | 基于互相关矩阵的扩展孔径二维doa估计方法 | |
CN110244258B (zh) | 用于双平行阵二维测向中扩展doa矩阵方法 | |
CN113673317B (zh) | 基于原子范数最小化可降维的二维离格doa估计方法 | |
CN109696657B (zh) | 一种基于矢量水听器的相干声源定位方法 | |
CN111965591A (zh) | 一种基于四阶累积量矢量化dft的测向估计方法 | |
CN108872930B (zh) | 扩展孔径二维联合对角化doa估计方法 | |
CN112255629A (zh) | 基于联合uca阵列的序贯esprit二维不相干分布源参数估计方法 | |
CN108267712A (zh) | 一种基于压缩平移互素阵列的doa估计方法及装置 | |
CN113552532B (zh) | 基于耦合张量分解的l型互质阵列波达方向估计方法 | |
CN113075610B (zh) | 一种基于互质极化阵列的差分阵列内插的doa估计方法 | |
WO2022166477A1 (zh) | 定位方法、装置、基站、计算机设备和存储介质 | |
CN114371441A (zh) | 虚拟阵列波达方向估计方法、装置、产品及存储介质 | |
CN113820655A (zh) | 一种基于Toeplitz矩阵重构和矩阵填充的互质阵相干信号DOA估计方法 | |
WO2021068494A1 (zh) | 基于平面互质阵列虚拟域张量空间谱搜索的高分辨精确二维波达方向估计方法 | |
CN114397619A (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 |