CN108594166A - 一种二维波达方向估计方法及装置 - Google Patents

一种二维波达方向估计方法及装置 Download PDF

Info

Publication number
CN108594166A
CN108594166A CN201810354351.0A CN201810354351A CN108594166A CN 108594166 A CN108594166 A CN 108594166A CN 201810354351 A CN201810354351 A CN 201810354351A CN 108594166 A CN108594166 A CN 108594166A
Authority
CN
China
Prior art keywords
matrix
angle parameter
dimension
noise
angle
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.)
Granted
Application number
CN201810354351.0A
Other languages
English (en)
Other versions
CN108594166B (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.)
Guangdong University of Technology
Original Assignee
Guangdong 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 Guangdong University of Technology filed Critical Guangdong University of Technology
Priority to CN201810354351.0A priority Critical patent/CN108594166B/zh
Publication of CN108594166A publication Critical patent/CN108594166A/zh
Application granted granted Critical
Publication of CN108594166B publication Critical patent/CN108594166B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

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/14—Systems for determining direction or deviation from predetermined direction

Landscapes

  • Physics & Mathematics (AREA)
  • Engineering & Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Remote Sensing (AREA)
  • Variable-Direction Aerials And Aerial Arrays (AREA)

Abstract

本发明公开了一种二维波达方向估计方法，包括：接收目标信号后，根据目标天线阵列和目标信号构造信号子空间对应的阵列流形，并基于阵列流形生成互协方差矩阵和传播算子，通过ESPRIT算法、互协方差矩阵和传播算子估计第一角度参数；同时构造噪声子空间对应的噪声投影矩阵，并基于噪声投影矩阵和阵列流形生成第一角度参数和第二角度参数的关系式；根据第一角度参数和关系式估计第二角度参数，进而估计波达方向。可见，所述方法基于第一角度参数估计第二角度参数，避免了二维角度参数的匹配，从而降低了计算的复杂度，同时具有良好的通用性。相应地，本发明公开的一种二维波达方向估计装置，也同样具有上述技术效果。

Description

一种二维波达方向估计方法及装置
技术领域
本发明涉及信号处理技术领域，更具体地说，涉及一种二维波达方向估计方法及装置。
背景技术
随着社会信息交流需求的急剧增加和个人移动通信的迅速普及，信号接收技术迎来了巨大挑战。
为了适应频谱效率和能量效率的数量级提升，大规模MIMO技术做为 MIMO技术的扩展和延伸技术已被广泛利用。其基本特征是在基站侧配置大规模的天线阵列(所包含的天线数量从几十至几千个)，显著提高频谱效率和能量效率。此外，在大规模MIMO系统中，使用简单的线性预编码和检测方法，噪声和快速衰落对系统的影响将逐渐消失，因此小区内干扰也得到了有效的抑制。这些优势使得大规模MIMO系统成为下一代5G移动通信系统的一大潜在关键技术。
在大规模MIMO系统中，波达方向估计成为重要研究方向之一。波达方向(Directionofarrival,DOA)是指空间信号的到达方向(各个信号到达阵列参考阵元的方向角，简称波达方向)。传统的高分辨率DOA估计算法由于需要对相关矩阵进行特征值分解得到噪声子空间，计算量巨大，难以应用于大规模MIMO系统中。
继而，传播子(Propagator Method)算法被提出，其利用线性运算代替特征分解，基于传播子的波达方向估计算法具有低复杂度特点，在大规模MIMO 的天线阵列及信号处理方面具有极大的潜力。YANG XIAO等提到一种低复杂度的二维双传播子算法(简称2P-PM算法)。该算法通过构造互相关矩阵来代替接收数据矩阵，从而抑制噪声分量带来的影响，并通过双传播子独立地估计二维角度参数，很大程度地降低算法计算量。但是该算法由于两个维度的角度参数由各自独立的ESPRIT步骤求出，故需要进一步的二维角度参数匹配，而匹配的计算量较大，一定程度上又增加了计算的复杂度；并且，在低信噪比下，2P-PM算法的估计性能大大降低，同时该算法要求信源个数小于某一维度上天线个数，这很大程度上限制了该算法的应用。
因此，如何降低波达方向估计的计算复杂度，提高估计性能，是本领域技术人员需要解决的问题。
发明内容
本发明的目的在于提供一种二维波达方向估计方法及装置，以实现降低波达方向估计的计算复杂度，提高估计性能。
为实现上述目的，本发明实施例提供了如下技术方案：
一种二维波达方向估计方法，包括：
通过目标天线阵列接收目标信号；
根据所述目标天线阵列和所述目标信号，构造信号子空间对应的阵列流形，并基于所述阵列流形生成互协方差矩阵和传播算子，通过ESPRIT算法、所述互协方差矩阵和所述传播算子估计第一角度参数；
根据所述目标天线阵列和所述目标信号，构造噪声子空间对应的噪声投影矩阵，并基于所述噪声投影矩阵和所述阵列流形生成所述第一角度参数和第二角度参数的关系式；
根据所述第一角度参数和所述关系式估计所述第二角度参数，并基于所述第一角度参数和所述第二角度参数估计波达方向。
其中，所述通过目标天线阵列接收目标信号，包括：
通过目标天线阵列接收目标信号X，X＝AS+N，其中，S为维度为(K×L) 的信源波形矩阵，N为维度为(MN×L)的高斯白噪声矩阵，A为所述阵列流形。
其中，所述根据所述目标天线阵列和所述目标信号，构造信号子空间对应的阵列流形，并基于所述阵列流形生成互协方差矩阵和传播算子，通过 ESPRIT算法、所述互协方差矩阵和所述传播算子估计第一角度参数，包括：
根据所述目标天线阵列和所述目标信号，建立沿X轴的方向矢量并基于所述方向矢量构造所述阵列流形A，其中， N为沿X轴等距分布的天线个数，M为沿Z轴等距分布的天线个数，Vi为所述第一角度参数，Ui为所述第二角度参数，K为信源个数；
将所述阵列流形A分为维度为((MN-n)×K)的矩阵Ap，其中，K＜n＜MN，和维度为(n×K)的矩阵Bp，并基于所述Ap和所述Bp分别对应的接收信号生成所述互协方差矩阵R21，R21＝[(BpS)(ApS)H]/L，其中，L为快拍数，S为所述维度为(K×L)的信源波形矩阵；
将所述Ap分为维度为(K×K)的矩阵Ap1和维度为((MN-n-K)×K)的矩阵 Ap2，并根据所述Ap的线性无关性构造所述Ap1和所述Ap2的线性关系式 Ap2＝PHAP1，其中，P为传播算子，维度为(K×(MN-n-K))；
将所述R21分为维度为(MN-n)×K)的矩阵R21a，和维度为((MN-n)×(n-K))的矩阵R21b，并基于所述R21a和R21b生成所述传播算子P，其中，
定义其中，IK为维度为K的单位矩阵，则其中，PC1和PC2分别为PC的前(MN-n-N)行和后(MN-n-N)行，Aa和Ab分别为AP的前(MN-n-N)行和后(MN-n-N)行，diag[]为对角矩阵，均为对角元素；
基于构造PC2AP1＝PC1AP1φZ、对进行特征值分解，并通过ui＝angle(φZ(i,i))λ/2πd估计所述第一角度参数 Ui。
其中，所述根据所述目标天线阵列和所述目标信号，构造噪声子空间对应的噪声投影矩阵，并基于所述噪声投影矩阵和所述阵列流形生成所述第一角度参数和第二角度参数的关系式，包括：
根据所述目标天线阵列和所述目标信号，构造所述噪声子空间对应的噪声矩阵Q和噪声投影矩阵Q0，令n＝MN/2，则Q0＝Q(QHQ)-1QH，其中，矩阵Q的维度为MN/2×(MN/2-K)，矩阵QHQ的维度为(MN/2-K)2，
由于P的维度为K×(MN/2-K)，PPH的维度为K×K，那么， E1＝PPH，E2＝(IK-PPH)-1，E＝IMN/2-K-PH(IK+PPH)-1P；
基于所述噪声投影矩阵Q0和所述阵列流形A构造关系式并基于所述构造所述第一角度参数和第二角度参数的关系式，其中，所述关系为：
其中，所述根据所述第一角度参数和所述关系式估计所述第二角度参数，并基于所述第一角度参数和所述第二角度参数估计波达方向，包括：
定义则表示为通过和所述第一角度参数估计所述第二角度参数，
基于信源个数构造目标多项式，利用所述目标多项式的最接近单位圆的根估计所述第二角度参数，并基于所述第一角度参数和所述第二角度参数估计波达方向。
一种二维波达方向估计装置，包括：
接收模块，用于通过目标天线阵列接收目标信号；
第一估计模块，用于根据所述目标天线阵列和所述目标信号，构造信号子空间对应的阵列流形，并基于所述阵列流形生成互协方差矩阵和传播算子，通过ESPRIT算法、所述互协方差矩阵和所述传播算子估计第一角度参数；
生成模块，用于根据所述目标天线阵列和所述目标信号，构造噪声子空间对应的噪声投影矩阵，并基于所述噪声投影矩阵和所述阵列流形生成所述第一角度参数和第二角度参数的关系式；
第二估计模块，用于根据所述第一角度参数和所述关系式估计所述第二角度参数，并基于所述第一角度参数和所述第二角度参数估计波达方向。
其中，所述接收模块具体用于：
通过目标天线阵列接收目标信号X，X＝AS+N，其中，S为维度为(K×L) 的信源波形矩阵，N为维度为(MN×L)的高斯白噪声矩阵，A为所述阵列流形。
其中，所述第一估计模块具体用于：
根据所述目标天线阵列和所述目标信号，建立沿X轴的方向矢量并基于所述方向矢量构造所述阵列流形A，其中， N为沿X轴等距分布的天线个数，M为沿Z轴等距分布的天线个数，Vi为所述第一角度参数，Ui为所述第二角度参数，K为信源个数；
将所述阵列流形A分为维度为((MN-n)×K)的矩阵Ap，其中，K＜n＜MN，和维度为(n×K)的矩阵Bp，并基于所述Ap和所述Bp分别对应的接收信号生成所述互协方差矩阵R21，R21＝[(BpS)(ApS)H]/L，其中，L为快拍数，S为所述维度为(K×L)的信源波形矩阵；
将所述Ap分为维度为(K×K)的矩阵Ap1和维度为((MN-n-K)×K)的矩阵 Ap2，并根据所述Ap的线性无关性构造所述Ap1和所述Ap2的线性关系式 Ap2＝PHAP1，其中，P为传播算子，维度为(K×(MN-n-K))；
将所述R21分为维度为(MN-n)×K)的矩阵R21a，和维度为 ((MN-n)×(n-K))的矩阵R21b，并基于所述R21a和R21b生成所述传播算子P，其中，
定义其中，IK为维度为K的单位矩阵，则其中，PC1和PC2分别为PC的前(MN-n-N)行和后(MN-n-N)行，Aa和Ab分别为AP的前(MN-n-N)行和后(MN-n-N)行，diag[]为对角矩阵，均为对角元素；
基于构造PC2AP1＝PC1AP1φZ、对进行特征值分解，并通过ui＝angle(φZ(i,i))λ/2πd估计所述第一角度参数 Ui。
其中，所述生成模块具体用于：
根据所述目标天线阵列和所述目标信号，构造所述噪声子空间对应的噪声矩阵Q和噪声投影矩阵Q0，令n＝MN/2，则Q0＝Q(QHQ)-1QH，其中，矩阵Q的维度为MN/2×(MN/2-K)，矩阵QHQ的维度为(MN/2-K)2，
由于P的维度为 K×(MN/2-K)，PPH的维度为K×K，那么， E1＝PPH，E2＝(IK-PPH)-1，E＝IMN/2-K-PH(IK+PPH)-1P；
基于所述噪声投影矩阵Q0和所述阵列流形A构造关系式并基于所述构造所述第一角度参数和第二角度参数的关系式，其中，所述关系为：
其中，所述第二估计模块具体用于：
定义则表示为通过和所述第一角度参数估计所述第二角度参数，
基于信源个数构造目标多项式，利用所述目标多项式的最接近单位圆的根估计所述第二角度参数，并基于所述第一角度参数和所述第二角度参数估计波达方向。
通过以上方案可知，本发明实施例提供的一种二维波达方向估计方法，包括：通过目标天线阵列接收目标信号；根据所述目标天线阵列和所述目标信号，构造信号子空间对应的阵列流形，并基于所述阵列流形生成互协方差矩阵和传播算子，通过ESPRIT算法、所述互协方差矩阵和所述传播算子估计第一角度参数；根据所述目标天线阵列和所述目标信号，构造噪声子空间对应的噪声投影矩阵，并基于所述噪声投影矩阵和所述阵列流形生成所述第一角度参数和第二角度参数的关系式；根据所述第一角度参数和所述关系式估计所述第二角度参数，并基于所述第一角度参数和所述第二角度参数估计波达方向。
可见，所述方法首先基于信号子空间构造阵列流形，并基于阵列流形生成互协方差矩阵和传播算子，通过ESPRIT算法、互协方差矩阵和传播算子估计出第一角度参数，进而构造噪声投影矩阵，利用噪声投影矩阵和阵列流形的正交性生成第一角度参数和第二角度参数的关系式，最后基于估计出的第一角度参数和所述关系式估计第二角度参数，并估计出波达方向。其中，所述互协方差矩阵能够有效地抑制噪声分量，提高估计性能；所述噪声投影矩阵可以降低计算的复杂度，同时基于第一角度参数估计第二角度参数，避免了二维角度参数匹配，从而降低了计算的复杂度，也避免了低信噪比下匹配错误的问题；并且，所述方法不受限于信源个数和某一维度上的天线个数，具有良好的通用性。
相应地，本发明实施例提供的一种二维波达方向估计装置，也同样具有上述技术效果。
附图说明
为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。
图1为本发明实施例公开的一种二维波达方向估计方法流程图；
图2为本发明实施例公开的(K<M或N时)PM Root MUSIC算法的二维角估计图；
图3为本发明实施例公开的(K>M或N时)PM Root MUSIC算法的二维角估计图；
图4为本发明实施例公开的PM Root MUSIC算法与2P-PM算法的均方差根对比图；
图5为本发明实施例公开的PM Root MUSIC算法与2P-PM算法的执行时间对比图；
图6为本发明实施例公开的一种二维波达方向估计装置示意图。
具体实施方式
下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。
本发明实施例公开了一种二维波达方向估计方法及装置，以实现降低波达方向估计的计算复杂度，提高估计性能。
参见图1，本发明实施例提供的一种二维波达方向估计方法，包括：
S101、通过目标天线阵列接收目标信号；
具体的，所述通过目标天线阵列接收目标信号，包括：通过目标天线阵列接收目标信号X，X＝AS+N，其中，S为维度为(K×L)的信源波形矩阵，N 为维度为(MN×L)的高斯白噪声矩阵，A为所述阵列流形。
S102、根据所述目标天线阵列和所述目标信号，构造信号子空间对应的阵列流形，并基于所述阵列流形生成互协方差矩阵和传播算子，通过ESPRIT 算法、所述互协方差矩阵和所述传播算子估计第一角度参数；
其中，所述根据所述目标天线阵列和所述目标信号，构造信号子空间对应的阵列流形，并基于所述阵列流形生成互协方差矩阵和传播算子，通过 ESPRIT算法、所述互协方差矩阵和所述传播算子估计第一角度参数，包括：
根据所述目标天线阵列和所述目标信号，建立沿X轴的方向矢量并基于所述方向矢量构造所述阵列流形A，其中， N为沿X轴等距分布的天线个数，M为沿Z轴等距分布的天线个数，Vi为所述第一角度参数，Ui为所述第二角度参数，K为信源个数；
将所述阵列流形A分为维度为((MN-n)×K)的矩阵Ap，其中，K＜n＜MN，和维度为(n×K)的矩阵Bp，并基于所述Ap和所述Bp分别对应的接收信号生成所述互协方差矩阵R21，R21＝[(BpS)(ApS)H]/L，其中，L为快拍数，S为所述维度为(K×L)的信源波形矩阵，Ap和Bp分别对应两个天线子阵列，由于阵间无共同阵元，故互协方差矩阵能够有效抑制噪声分量；
将所述Ap分为维度为(K×K)的矩阵Ap1和维度为((MN-n-K)×K)的矩阵 Ap2，并根据所述Ap的线性无关性构造所述Ap1和所述Ap2的线性关系式 Ap2＝PHAP1，其中，P为传播算子，维度为(K×(MN-n-K))；
需要说明的是，假设信源为非相干信源，那么Ap的秩为K，且Ap的前K 行是线性无关的，因此可得到Ap1和Ap2的线性关系式Ap2＝PHAP1；
将所述R21分为维度为(MN-n)×K)的矩阵R21a，和维度为((MN-n)×(n-K))的矩阵R21b，并基于所述R21a和R21b生成所述传播算子P，其中，
需要说明的是，将公式Ap2＝PHAP1代入公式R21＝[(BpS)(ApS)H]/L，即可得到而由于()+为伪逆矩阵，故传播算子
为了避免谱峰搜索，采用ESPRIT(旋转不变因子)算法估计第一角度参数：
定义其中，IK为维度为K的单位矩阵，则其中，PC1和PC2分别为PC的前(MN-n-N)行和后(MN-n-N)行，Aa和Ab分别为AP的前(MN-n-N)行和后(MN-n-N)行，diag[]为对角矩阵，均为对角元素；
基于构造PC2AP1＝PC1AP1φZ、对进行特征值分解，并通过ui＝angle(φZ(i,i))λ/2πd估计所述第一角度参数 Ui；
需要说明的是，由于和φZ有相同的特征值，故可通过对进行特征值分解来估计第一角度参数Ui，在公式ui＝angle(φZ(i,i))λ/2πd 中，angle()为取括号内元素的相位。
S103、根据所述目标天线阵列和所述目标信号，构造噪声子空间对应的噪声投影矩阵，并基于所述噪声投影矩阵和所述阵列流形生成所述第一角度参数和第二角度参数的关系式；
其中，所述根据所述目标天线阵列和所述目标信号，构造噪声子空间对应的噪声投影矩阵，并基于所述噪声投影矩阵和所述阵列流形生成所述第一角度参数和第二角度参数的关系式，包括：
根据所述目标天线阵列和所述目标信号，构造所述噪声子空间对应的噪声矩阵Q和噪声投影矩阵Q0，令n＝MN/2，则Q0＝Q(QHQ)-1QH，其中，矩阵Q的维度为MN/2×(MN/2-K)，矩阵QHQ的维度为(MN/2-K)2，其中，令n＝MN/2仅是本方法中n值的一个取值代表，故n的值可根据实际情况进行变换；噪声投影矩阵Q0可以降低计算的复杂度，以提高估计性能；
由于P的维度为 K×(MN/2-K)，PPH的维度为K×K，那么， E1＝PPH，E2＝(IK-PPH)-1，E＝IMN/2-K-PH(IK+PPH)-1P；
基于所述噪声投影矩阵Q0和所述阵列流形A构造关系式并基于所述构造所述第一角度参数和第二角度参数的关系式，其中，所述关系为：
S104、根据所述第一角度参数和所述关系式估计所述第二角度参数，并基于所述第一角度参数和所述第二角度参数估计波达方向。
其中，所述根据所述第一角度参数和所述关系式估计所述第二角度参数，并基于所述第一角度参数和所述第二角度参数估计波达方向，包括：
定义则表示为通过和所述第一角度参数估计所述第二角度参数，
基于信源个数构造目标多项式，利用所述目标多项式的最接近单位圆的根估计所述第二角度参数，并基于所述第一角度参数和所述第二角度参数估计波达方向。
具体的，可构造多项式：其中，关于Z的2(N-1)阶多项式可以表示为：
其中，Ci表示Zi的系数，di对应Di每条对角线上元素之和，定义矩阵Di中第x行，y列上的元素为dxy，可知
例如：当非相干信源个数为3时，即N＝3，则有：
由此可得出下述多项式方程：
求解所述方程，取最接近单位圆的根，进而估计第二角度参数，其中，由于第二角度参数是基于第一角度参数而得，故第一角度参数和第二角度参数具有意义对应性，即实现了角度参数的匹配。
可见，本实施例提供的一种二维波达方向估计方法首先基于信号子空间构造阵列流形，并基于阵列流形生成互协方差矩阵和传播算子，通过ESPRIT 算法、互协方差矩阵和传播算子估计出第一角度参数，进而构造噪声投影矩阵，利用噪声投影矩阵和阵列流形的正交性生成第一角度参数和第二角度参数的关系式，最后基于估计出的第一角度参数和所述关系式估计第二角度参数，并估计出波达方向。其中，所述互协方差矩阵能够有效地抑制噪声分量，提高估计性能；所述噪声投影矩阵可以降低计算的复杂度，同时基于第一角度参数估计第二角度参数，避免了二维角度参数匹配，从而降低了计算的复杂度，也避免了低信噪比下匹配错误的问题；并且，所述方法不受限于信源个数和某一维度上的天线个数，具有良好的通用性。
基于上述实施例，需要说明的是，所述基于所述噪声投影矩阵和所述阵列流形生成所述第一角度参数和第二角度参数的关系式的计算过程可优化。
由于矩阵Q的维度为MN/2×(MN/2-K)，矩阵QHQ的维度为(MN/2-K)2, 那么(QHQ)-1的计算复杂度为：O(2(MN/2-K)2(MN/2)+(MN/2-K)3)，进而可知： Q0＝Q(QHQ)-1QH的计算复杂度为：
O(2(MN/2-K)2(MN/2)+(MN/2-K)3+(MN/2)2(MN/2-K))；
对上述计算过程进行优化，则：
由于P的维度为K×(MN/2-K)，PPH的维度为K×K，则：
IMN/2-K-PH(PPH+IK)-1P的计算复杂度为：O(3K2(MN/2-K)+K3)；
而故优化后的计算复杂度为：
O(3K2(MN/2-K)+5K3+(MN/2-K)2K)，而由于K远小于MN，可见优化后的计算复杂度大大降低。
将阵列流形A的子矩阵AP表示为：
其中，为克罗内克积；由此可得到第一角度参数和第二角度参数的关系式：
具体的，由于可知矩阵QH的列所形成的信号空间与阵列流形子矩阵Ap的列所形成的信号空间正交，即阵列流形Ap任意方向向量d(θ)与噪声子空间Q的列相互正交，即QHAP的模的平方为0，其中，QH为Q的共轭转置矩阵。
需要说明的是，在实际理论中，噪声矩阵Q与Ap正交，由于噪声存在， Q与Ap非严格正交，故通过投影变换得到噪声投影矩阵Q0，使其与Ap严格正交，实际计算过程中则利用Q0的共轭转置矩阵Q0 H进行相应计算，即可构造公式至此，便可利用Q0与Ap正交性生成第一角度参数和第二角度参数的关系式。其中，噪声矩阵Q经过投影变换可得到噪声投影矩阵Q0。
基于上述实施例，需要说明的是，所述角度参数为仰角和方位角，其求解顺序可互换。即：可先求出仰角，再依据仰角和关系式求解方位角；也可先求出方位角，再依据方位角和关系式求解仰角。
基于上述实施例，需要说明的是，本说明书提供的二维波达方向估计方法结合了2P-PM算法和Root MUSIC算法，故本说明书提供的二维波达方向估计方法可称为PM RootMUSIC算法。
基于上述实施例，利用仿真实验验证上述二维波达方向估计方法。需要说明的是，可通过修改n值来构造不同的互协方差矩阵。
本仿真实验令n＝MN/2，通过仿真实验，比较在远场下本发明提供的方法与2P-PM算法的波达方向估计性能以及算法运行时间。
仿真初始参数如下，仿真实验次数为500次，相邻天线间距为半波长，矩形天线阵列规模为M＝N＝16(即16*16根天线)，快拍数L为100。考虑三个窄带信号作为信源，其仰角和方位角分别为(80，70)，(70，50)，(60， 30)。定义均方根误差(RMSE)为：
其中，K为信号源个数，S为蒙特卡洛实验次数.和分别为仰角和方位角的估计值。
参见图2，图2所示内容为信源个数K＝3，信噪比SNR＝15dB，天线阵列规模M＝N＝16情况下，本发明算法的二维角估计图，即当信源个数小于某一维度上的天线数量(K<M或N)时，本发明算法有效。
参见图3，图3所示内容为信源个数K＝9，信噪比SNR＝15dB，天线阵列规模M＝N＝8情况下，本发明算法的二维角估计图，即当信源个数大于某一维度上的天线数量(K>M或N)时，本发明算法仍然有效，即克服了信源个数和某一维度上的天线数量的限制；而2P-PM算法由于通过其中一列天线作为子阵列来估计方位角，对应的传播子P的维度为(K×(N-K))，当K>N时，2P-PM 算法失效。
参见图4，图4所示内容为信源个数K＝3，天线数量M＝N＝16，快拍数 L＝100，不同信噪比下，两种估计算法的估计性能。由此可见，本发明提供的方法的估计性能有明显的提高。此外，虽然2P-PM算法通过构建协方差来抑制噪声分量，但由于需要进行角度参数匹配，在低信噪比下会存在匹配错误等问题，造成在信噪比小于10dB时，其算法效果明显恶化。而本发明通过对噪声子空间进行投影变换以及自动匹配的方法，解决了在低信噪下，算法效果恶化问题。
参见图5，图5所示内容为信源个数K＝3，信噪比SNR＝10dB，快拍数 L＝100，不同阵列规模情况下，两种估计算法的执行时间。由此可见，本发明在天线规模为100-900范围内的执行时间与2P-PM算法接近，具备高效低复杂度的优点。
下面对本发明实施例提供的一种二维波达方向估计装置进行介绍，下文描述的一种二维波达方向估计装置与上文描述的一种二维波达方向估计方法可以相互参照。
参见图6，本发明实施例提供的一种二维波达方向估计装置，包括：
接收模块601，用于通过目标天线阵列接收目标信号；
第一估计模块602，用于根据所述目标天线阵列和所述目标信号，构造信号子空间对应的阵列流形，并基于所述阵列流形生成互协方差矩阵和传播算子，通过ESPRIT算法、所述互协方差矩阵和所述传播算子估计第一角度参数；
生成模块603，用于根据所述目标天线阵列和所述目标信号，构造噪声子空间对应的噪声投影矩阵，并基于所述噪声投影矩阵和所述阵列流形生成所述第一角度参数和第二角度参数的关系式；
第二估计模块604，用于根据所述第一角度参数和所述关系式估计所述第二角度参数，并基于所述第一角度参数和所述第二角度参数估计波达方向。
其中，所述接收模块具体用于：
通过目标天线阵列接收目标信号X，X＝AS+N，其中，S为维度为(K×L) 的信源波形矩阵，N为维度为(MN×L)的高斯白噪声矩阵，A为所述阵列流形。
其中，所述第一估计模块具体用于：
根据所述目标天线阵列和所述目标信号，建立沿X轴的方向矢量并基于所述方向矢量构造所述阵列流形A，其中， N为沿X轴等距分布的天线个数，M为沿Z轴等距分布的天线个数，Vi为所述第一角度参数，Ui为所述第二角度参数，K为信源个数；
将所述阵列流形A分为维度为((MN-n)×K)的矩阵Ap，其中，K＜n＜MN，和维度为(n×K)的矩阵Bp，并基于所述Ap和所述Bp分别对应的接收信号生成所述互协方差矩阵R21，R21＝[(BpS)(ApS)H]/L，其中，L为快拍数，S为所述维度为(K×L)的信源波形矩阵；
将所述Ap分为维度为(K×K)的矩阵Ap1和维度为((MN-n-K)×K)的矩阵 Ap2，并根据所述Ap的线性无关性构造所述Ap1和所述Ap2的线性关系式 Ap2＝PHAP1，其中，P为传播算子，维度为(K×(MN-n-K))；
将所述R21分为维度为(MN-n)×K)的矩阵R21a，和维度为((MN-n)×(n-K))的矩阵R21b，并基于所述R21a和R21b生成所述传播算子P，其中，
定义其中，IK为维度为K的单位矩阵，则其中，PC1和PC2分别为PC的前(MN-n-N)行和后(MN-n-N)行，Aa和Ab分别为AP的前(MN-n-N)行和后(MN-n-N)行，diag[]为对角矩阵，均为对角元素；
基于构造PC2AP1＝PC1AP1φZ、对进行特征值分解，并通过ui＝angle(φZ(i,i))λ/2πd估计所述第一角度参数 Ui。
其中，所述生成模块具体用于：
根据所述目标天线阵列和所述目标信号，构造所述噪声子空间对应的噪声矩阵Q和噪声投影矩阵Q0，令n＝MN/2，则Q0＝Q(QHQ)-1QH，其中，矩阵Q的维度为MN/2×(MN/2-K)，矩阵QHQ的维度为(MN/2-K)2，
由于P的维度为 K×(MN/2-K)，PPH的维度为K×K，那么， E1＝PPH，E2＝(IK-PPH)-1，E＝IMN/2-K-PH(IK+PPH)-1P；
基于所述噪声投影矩阵Q0和所述阵列流形A构造关系式并基于所述构造所述第一角度参数和第二角度参数的关系式，其中，所述关系为：
其中，所述第二估计模块具体用于：
定义则表示为通过和所述第一角度参数估计所述第二角度参数，
基于信源个数构造目标多项式，利用所述目标多项式的最接近单位圆的根估计所述第二角度参数，并基于所述第一角度参数和所述第二角度参数估计波达方向。
本说明书中各个实施例采用递进的方式描述，每个实施例重点说明的都是与其他实施例的不同之处，各个实施例之间相同相似部分互相参见即可。
对所公开的实施例的上述说明，使本领域专业技术人员能够实现或使用本发明。对这些实施例的多种修改对本领域的专业技术人员来说将是显而易见的，本文中所定义的一般原理可以在不脱离本发明的精神或范围的情况下，在其它实施例中实现。因此，本发明将不会被限制于本文所示的这些实施例，而是要符合与本文所公开的原理和新颖特点相一致的最宽的范围。

Claims (10)

1.一种二维波达方向估计方法，其特征在于，包括：
通过目标天线阵列接收目标信号；
根据所述目标天线阵列和所述目标信号，构造信号子空间对应的阵列流形，并基于所述阵列流形生成互协方差矩阵和传播算子，通过ESPRIT算法、所述互协方差矩阵和所述传播算子估计第一角度参数；
根据所述目标天线阵列和所述目标信号，构造噪声子空间对应的噪声投影矩阵，并基于所述噪声投影矩阵和所述阵列流形生成所述第一角度参数和第二角度参数的关系式；
根据所述第一角度参数和所述关系式估计所述第二角度参数，并基于所述第一角度参数和所述第二角度参数估计波达方向。
2.根据权利要求1所述的二维波达方向估计方法，其特征在于，所述通过目标天线阵列接收目标信号，包括：
通过目标天线阵列接收目标信号X，X＝AS+N，其中，S为维度为(K×L)的信源波形矩阵，N为维度为(MN×L)的高斯白噪声矩阵，A为所述阵列流形。
3.根据权利要求1所述的二维波达方向估计方法，其特征在于，所述根据所述目标天线阵列和所述目标信号，构造信号子空间对应的阵列流形，并基于所述阵列流形生成互协方差矩阵和传播算子，通过ESPRIT算法、所述互协方差矩阵和所述传播算子估计第一角度参数，包括：
根据所述目标天线阵列和所述目标信号，建立沿X轴的方向矢量并基于所述方向矢量构造所述阵列流形A，其中，N为沿X轴等距分布的天线个数，M为沿Z轴等距分布的天线个数，Vi为所述第一角度参数，Ui为所述第二角度参数，K为信源个数；
将所述阵列流形A分为维度为((MN-n)×K)的矩阵Ap，其中，K＜n＜MN，和维度为(n×K)的矩阵Bp，并基于所述Ap和所述Bp分别对应的接收信号生成所述互协方差矩阵R21，R21＝[(BpS)(ApS)H]/L，其中，L为快拍数，S为所述维度为(K×L)的信源波形矩阵；
将所述Ap分为维度为(K×K)的矩阵Ap1和维度为((MN-n-K)×K)的矩阵Ap2，并根据所述Ap的线性无关性构造所述Ap1和所述Ap2的线性关系式Ap2＝PHAP1，其中，P为传播算子，维度为(K×(MN-n-K))；
将所述R21分为维度为(MN-n)×K)的矩阵R21a，和维度为((MN-n)×(n-K))的矩阵R21b，并基于所述R21a和R21b生成所述传播算子P，其中，
定义其中，IK为维度为K的单位矩阵，则其中，PC1和PC2分别为PC的前(MN-n-N)行和后(MN-n-N)行，Aa和Ab分别为AP的前(MN-n-N)行和后(MN-n-N)行，diag[]为对角矩阵，均为对角元素；
基于构造PC2AP1＝PC1AP1φZ、对进行特征值分解，并通过ui＝angle(φZ(i,i))λ/2πd估计所述第一角度参数Ui。
4.根据权利要求1所述的二维波达方向估计方法，其特征在于，所述根据所述目标天线阵列和所述目标信号，构造噪声子空间对应的噪声投影矩阵，并基于所述噪声投影矩阵和所述阵列流形生成所述第一角度参数和第二角度参数的关系式，包括：
根据所述目标天线阵列和所述目标信号，构造所述噪声子空间对应的噪声矩阵Q和噪声投影矩阵Q0，令n＝MN/2，则Q0＝Q(QHQ)-1QH，其中，矩阵Q的维度为MN/2×(MN/2-K)，矩阵QHQ的维度为(MN/2-K)2，
由于P的维度为K×(MN/2-K)，PPH的维度为K×K，那么，E1＝PPH，E2＝(IK-PPH)-1，E＝IMN/2-K-PH(IK+PPH)-1P；
基于所述噪声投影矩阵Q0和所述阵列流形A构造关系式并基于所述构造所述第一角度参数和第二角度参数的关系式，其中，所述关系为：
5.根据权利要求1所述的二维波达方向估计方法，其特征在于，所述根据所述第一角度参数和所述关系式估计所述第二角度参数，并基于所述第一角度参数和所述第二角度参数估计波达方向，包括：
定义则表示为
通过和所述第一角度参数估计所述第二角度参数，
基于信源个数构造目标多项式，利用所述目标多项式的最接近单位圆的根估计所述第二角度参数，并基于所述第一角度参数和所述第二角度参数估计波达方向。
6.一种二维波达方向估计装置，其特征在于，包括：
接收模块，用于通过目标天线阵列接收目标信号；
第一估计模块，用于根据所述目标天线阵列和所述目标信号，构造信号子空间对应的阵列流形，并基于所述阵列流形生成互协方差矩阵和传播算子，通过ESPRIT算法、所述互协方差矩阵和所述传播算子估计第一角度参数；
生成模块，用于根据所述目标天线阵列和所述目标信号，构造噪声子空间对应的噪声投影矩阵，并基于所述噪声投影矩阵和所述阵列流形生成所述第一角度参数和第二角度参数的关系式；
第二估计模块，用于根据所述第一角度参数和所述关系式估计所述第二角度参数，并基于所述第一角度参数和所述第二角度参数估计波达方向。
7.根据权利要求6所述的二维波达方向估计装置，其特征在于，所述接收模块具体用于：
通过目标天线阵列接收目标信号X，X＝AS+N，其中，S为维度为(K×L)的信源波形矩阵，N为维度为(MN×L)的高斯白噪声矩阵，A为所述阵列流形。
8.根据权利要求6所述的二维波达方向估计装置，其特征在于，所述第一估计模块具体用于：
根据所述目标天线阵列和所述目标信号，建立沿X轴的方向矢量并基于所述方向矢量构造所述阵列流形A，其中，N为沿X轴等距分布的天线个数，M为沿Z轴等距分布的天线个数，Vi为所述第一角度参数，Ui为所述第二角度参数，K为信源个数；
将所述阵列流形A分为维度为((MN-n)×K)的矩阵Ap，其中，K＜n＜MN，和维度为(n×K)的矩阵Bp，并基于所述Ap和所述Bp分别对应的接收信号生成所述互协方差矩阵R21，R21＝[(BpS)(ApS)H]/L，其中，L为快拍数，S为所述维度为(K×L)的信源波形矩阵；
将所述Ap分为维度为(K×K)的矩阵Ap1和维度为((MN-n-K)×K)的矩阵Ap2，并根据所述Ap的线性无关性构造所述Ap1和所述Ap2的线性关系式Ap2＝PHAP1，其中，P为传播算子，维度为(K×(MN-n-K))；
将所述R21分为维度为(MN-n)×K)的矩阵R21a，和维度为((MN-n)×(n-K))的矩阵R21b，并基于所述R21a和R21b生成所述传播算子P，其中，
定义其中，IK为维度为K的单位矩阵，则其中，PC1和PC2分别为PC的前(MN-n-N)行和后(MN-n-N)行，Aa和Ab分别为AP的前(MN-n-N)行和后(MN-n-N)行，diag[]为对角矩阵，均为对角元素；
基于构造PC2AP1＝PC1AP1φZ、对进行特征值分解，并通过ui＝angle(φZ(i,i))λ/2πd估计所述第一角度参数Ui。
9.根据权利要求6所述的二维波达方向估计装置，其特征在于，所述生成模块具体用于：
根据所述目标天线阵列和所述目标信号，构造所述噪声子空间对应的噪声矩阵Q和噪声投影矩阵Q0，令n＝MN/2，则Q0＝Q(QHQ)-1QH，其中，矩阵Q的维度为MN/2×(MN/2-K)，矩阵QHQ的维度为(MN/2-K)2，
由于P的维度为K×(MN/2-K)，PPH的维度为K×K，那么，E1＝PPH，E2＝(IK-PPH)-1，E＝IMN/2-K-PH(IK+PPH)-1P；
基于所述噪声投影矩阵Q0和所述阵列流形A构造关系式并基于所述构造所述第一角度参数和第二角度参数的关系式，其中，所述关系为：
10.根据权利要求6所述的二维波达方向估计装置，其特征在于，所述第二估计模块具体用于：
定义则表示为
通过和所述第一角度参数估计所述第二角度参数，
基于信源个数构造目标多项式，利用所述目标多项式的最接近单位圆的根估计所述第二角度参数，并基于所述第一角度参数和所述第二角度参数估计波达方向。
CN201810354351.0A 2018-04-19 2018-04-19 一种二维波达方向估计方法及装置 Active CN108594166B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201810354351.0A CN108594166B (zh) 2018-04-19 2018-04-19 一种二维波达方向估计方法及装置

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201810354351.0A CN108594166B (zh) 2018-04-19 2018-04-19 一种二维波达方向估计方法及装置

Publications (2)

Publication Number Publication Date
CN108594166A true CN108594166A (zh) 2018-09-28
CN108594166B CN108594166B (zh) 2022-03-25

Family

ID=63613931

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201810354351.0A Active CN108594166B (zh) 2018-04-19 2018-04-19 一种二维波达方向估计方法及装置

Country Status (1)

Country Link
CN (1) CN108594166B (zh)

Cited By (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN109709510A (zh) * 2018-12-24 2019-05-03 贵州航天计量测试技术研究所 一种相干信源二维波达方向的估计方法及系统
CN110109053A (zh) * 2019-04-02 2019-08-09 华南理工大学 一种未知声速环境下快速doa估计方法
CN110907926A (zh) * 2019-10-29 2020-03-24 长江大学 基于传播算子的双基地emvs-mimo雷达快速目标定位算法及装置
CN114741649A (zh) * 2022-03-28 2022-07-12 电子科技大学 一种基于sdae-dnn的低信噪比波达方向估计方法
CN116054909A (zh) * 2022-12-16 2023-05-02 广东工业大学 一种基于改进型传播算子的卫星通信抗干扰方法

Citations (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2009192516A (ja) * 2008-01-16 2009-08-27 Mitsubishi Electric Corp 伝搬遅延時間測定装置及びレーダ装置
CN102169170A (zh) * 2010-12-29 2011-08-31 电子科技大学 一种相干分布式信号二维波达角的测定方法
CN104459740A (zh) * 2014-11-12 2015-03-25 广东工业大学 一种定位终端的高精度位置差分定位方法
CN106500698A (zh) * 2016-10-31 2017-03-15 南京工程学院 一种基于传播算子的室内定位方法及装置
CN106526530A (zh) * 2016-09-30 2017-03-22 天津大学 基于传播算子的2‑l型阵列二维doa估计算法
CN106707228A (zh) * 2016-11-30 2017-05-24 广东工业大学 一种结合四阶累积量与盖尔圆改进的信号源个数估计方法
CN106842115A (zh) * 2017-01-25 2017-06-13 天津大学 基于root‑music算法及秩损原理的二维测向方法
CN107884741A (zh) * 2017-10-30 2018-04-06 北京理工大学 一种多球阵列多宽带声源快速定向方法

Patent Citations (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2009192516A (ja) * 2008-01-16 2009-08-27 Mitsubishi Electric Corp 伝搬遅延時間測定装置及びレーダ装置
CN102169170A (zh) * 2010-12-29 2011-08-31 电子科技大学 一种相干分布式信号二维波达角的测定方法
CN104459740A (zh) * 2014-11-12 2015-03-25 广东工业大学 一种定位终端的高精度位置差分定位方法
CN106526530A (zh) * 2016-09-30 2017-03-22 天津大学 基于传播算子的2‑l型阵列二维doa估计算法
CN106500698A (zh) * 2016-10-31 2017-03-15 南京工程学院 一种基于传播算子的室内定位方法及装置
CN106707228A (zh) * 2016-11-30 2017-05-24 广东工业大学 一种结合四阶累积量与盖尔圆改进的信号源个数估计方法
CN106842115A (zh) * 2017-01-25 2017-06-13 天津大学 基于root‑music算法及秩损原理的二维测向方法
CN107884741A (zh) * 2017-10-30 2018-04-06 北京理工大学 一种多球阵列多宽带声源快速定向方法

Non-Patent Citations (4)

* Cited by examiner, † Cited by third party
Title
XIAO YANG,ET AL: "A New Low Complexity DOA estimation algorithm for Massive MIMO Systems", 《2016 IEEE INTERNATIONAL CONFERENCE ON CONSUMER ELECTRONICS-CHINA (ICCE-CHINA)》 *
YANG-YANG DONG，ET AL: "Computationally Efficient 2-D DOA Estimation for L-Shaped Array With Automatic Pairing", 《IEEE ANTENNAS AND WIRELESS PROPAGATION LETTERS》 *
李小宇: "大规模MIMO相关技术研究", 《中国优秀硕士学位论文全文数据库 信息科技辑》 *
杨骁: "大规模MIMO系统中低复杂度DOA估计算法的研究", 《中国优秀硕士学位论文全文数据库 信息科技辑》 *

Cited By (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN109709510A (zh) * 2018-12-24 2019-05-03 贵州航天计量测试技术研究所 一种相干信源二维波达方向的估计方法及系统
CN110109053A (zh) * 2019-04-02 2019-08-09 华南理工大学 一种未知声速环境下快速doa估计方法
CN110907926A (zh) * 2019-10-29 2020-03-24 长江大学 基于传播算子的双基地emvs-mimo雷达快速目标定位算法及装置
CN114741649A (zh) * 2022-03-28 2022-07-12 电子科技大学 一种基于sdae-dnn的低信噪比波达方向估计方法
CN116054909A (zh) * 2022-12-16 2023-05-02 广东工业大学 一种基于改进型传播算子的卫星通信抗干扰方法
CN116054909B (zh) * 2022-12-16 2023-08-08 广东工业大学 一种基于改进型传播算子的卫星通信抗干扰方法

Also Published As

Publication number Publication date
CN108594166B (zh) 2022-03-25

Similar Documents

Publication Publication Date Title
CN108594166A (zh) 一种二维波达方向估计方法及装置
CN106972882B (zh) 基于虚拟域空间功率谱估计的互质阵列自适应波束成形方法
CN104977558B (zh) 一种基于贝叶斯压缩感知的分布源中心波达方向估计方法
CN113472705B (zh) 基于zc序列的ris辅助信道的估计及预测方法
CN107589399B (zh) 基于多采样虚拟信号奇异值分解的互质阵列波达方向估计方法
CN104991236B (zh) 一种单基地mimo雷达非圆信号相干源波达方向估计方法
CN107703478B (zh) 基于互相关矩阵的扩展孔径二维doa估计方法
CN107561484A (zh) 基于内插互质阵列协方差矩阵重建的波达方向估计方法
CN104698430A (zh) 一种用于提高精度的基于虚拟天线阵列的角度估计方法
CN104535987A (zh) 适用于均匀圆阵列声纳系统的幅相误差自校正方法
CN108872930B (zh) 扩展孔径二维联合对角化doa估计方法
CN111352063A (zh) 一种均匀面阵中基于多项式求根的二维测向估计方法
CN101483280B (zh) 一种稳健波束合成器权值求解方法
CN112255629A (zh) 基于联合uca阵列的序贯esprit二维不相干分布源参数估计方法
CN103323810B (zh) 一种l阵方位角和俯仰角配对的信号处理方法
CN113835063A (zh) 一种无人机阵列幅相误差与信号doa联合估计方法
CN105306117A (zh) 一种基于协方差矩阵扩展的半虚拟天线阵波束形成方法
Xu et al. Tensor‐Based Angle and Range Estimation Method in Monostatic FDA‐MIMO Radar
CN104459627B (zh) 基于联合交替优化的降秩波束形成方法
CN107479025B (zh) 一种单快拍大规模线阵空间谱估计方法
CN114520699B (zh) 可重构智能表面辅助毫米波通信的低复杂度信道估计方法
CN111614386B (zh) 一种快速搜索方法
CN101557250B (zh) 一种基于方向预测的智能天线权值生成方法和装置
Cong et al. Feasible sparse spectrum fitting of DOA and range estimation for collocated FDA-MIMO radars
CN113420411B (zh) 用于无线信号的高分辨率窄带doa估计算法及实现方法

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