CN111722178A - 一种基于指向性模型数值求解的远场窄带信号来波方向估计方法 - Google Patents

一种基于指向性模型数值求解的远场窄带信号来波方向估计方法 Download PDF

Info

Publication number
CN111722178A
CN111722178A CN202010474190.6A CN202010474190A CN111722178A CN 111722178 A CN111722178 A CN 111722178A CN 202010474190 A CN202010474190 A CN 202010474190A CN 111722178 A CN111722178 A CN 111722178A
Authority
CN
China
Prior art keywords
far
incoming wave
wave direction
signal
band signal
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
CN202010474190.6A
Other languages
English (en)
Other versions
CN111722178B (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.)
Southeast University
Original Assignee
Southeast University
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 Southeast University filed Critical Southeast University
Priority to CN202010474190.6A priority Critical patent/CN111722178B/zh
Publication of CN111722178A publication Critical patent/CN111722178A/zh
Application granted granted Critical
Publication of CN111722178B publication Critical patent/CN111722178B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

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
    • G01S3/143—Systems for determining direction or deviation from predetermined direction by vectorial combination of signals derived from differently oriented antennae
    • 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/78—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 electromagnetic waves other than radio waves
    • 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/80—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 ultrasonic, sonic or infrasonic waves
    • 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)
  • Electromagnetism (AREA)
  • Measurement Of Velocity Or Position Using Acoustic Or Ultrasonic Waves (AREA)

Abstract

本发明公开了一种基于指向性模型数值求解的远场窄带信号来波方向估计方法，该方法包括：第一步：获取M通道阵元数据x1(n),x2(n),…,xM(n)；第二步：根据给定的阵列参数自动生成等余弦间隔的预引导波束扫描角θk和总的预引导波束个数K；第三步：对已知频率f0的远场窄带信号进行频域波束形成，并求得该已知频率远场窄带信号的各预引导波束扫描角所对应的波束功率谱P(θk)；第四步：搜索波束功率谱P(θk)最大值所对应的预引导波束号kp；第五步：利用二分数值求解法求解信号来波方向相对偏差
Figure DDA0002515294780000011
第六步，利用
Figure DDA0002515294780000012
估计远场窄带信号的来波方向。该算法基于常规波束形成，实现简单，运算量小，且估计精度高，适用于实时的工程应用场合。

Description

一种基于指向性模型数值求解的远场窄带信号来波方向估计 方法
技术领域
本发明属于信号处理领域，尤其涉及一种基于指向性模型数值求解的远场窄带信号来波方向估计方法。
背景技术
目标来波方向估计是雷达、声纳、声学、语音和无线通信等阵列信号处理应用中的一个基本问题，尤其在水声和电子侦察处理中扮演了极其重要的角色。此外，在现代移动通信系统中，来波方向估计在提高无线系统性方面有着重要应用。
目前国内外学者提出了很多目标来波方向估计方法，如最小方差无失真响应、基于特征向量的来波方向估计算法(如MUSIC、Root-MUSIC、RV-Root-MUSIC)。这些算法通常能达到很高的精度。然而，这些基于峰值搜索的算法需要很大的计算复杂度。此外还存在着估计性能低、复杂度高和相位模糊等问题。同时，这些算法对信号失配非常敏感。
在波束形成方法中，常规波束形成是最常用的一种非自适应波束形成方法。常规波束形成使用一组固定的加权和时间延迟对来自阵列中传感器的信号进行组合叠加，这个叠加过程主要是使用关于传感器在空间中的位置和感兴趣的信号方向的信息。常规波束形成的优点是它对假定的和实际的信号波前之间的信号失配具有较强的鲁棒性，并且只需要少量样本即可稳健地估计目标来波方向，易于实现，因此在单目标场景中有着较好的性能。目前基于常规波束形成的来波方向估计方法主要有：峰值搜索算法和抛物线插值算法。
谱峰搜索法利用常规波束形成结果使角谱最大化，并返回M个最大值的位置作为M个目标的来波方向估计。谱峰搜索法对于白噪声中单个远场、窄带和电源的来波方向估计而言，是已知的最佳地最大似然估计方法。谱峰搜索法的来波方向估计精度取决于预引导波束扫描角的扫描间隔。为了获得更高的估计精度，需要更多的预引导波束扫描角。然而，更多的预引导波束扫描角将带来更多的计算负担。此外，受方位估计瑞利极限的限制，仅通过增加预引导波束扫描角的个数是无法无限提高来波方向估计精度的。
抛物线插值算法是基于对主瓣角谱的抛物线拟合，易于实现。对于大多数目标的来波方向估计，估计精度较高。然而，对于某些来波方向，角谱的主瓣与抛物线模型不匹配，在这种情况下，来波方向估计精度下降。
发明内容
技术问题：本发明提出了一种利用给定均匀线阵的理论波束图的模型求解算法。该算法是基于二分法的数值解，采用二分法数值求解信号来波方向相对偏差从而精确地估计远场窄带信号的来波方向。由于采用的是理论波束模式模型，因此不存在模型失配。与抛物线插值算法相比，该算法对信号来波方向相对偏差具有更强的鲁棒性，并且具有更好的估计性能。
技术方案：本发明的一种基于指向性模型数值求解远场窄带信号来波方向估计方法，包括以下步骤：
第一步：获取待处理的均匀线阵列M通道阵元远场窄带信号采样数据序列x1(n)，x2(n),…,xM(n),n＝0,1,…,N-1：从所述阵元基元接收N个采样点的实时采集数据作为待处理的数据序列x1(n),x2(n),…,xM(n),n＝0,1,…,N-1；或从存储器中提取从某一时刻开始接收到的N个信号采样点的数据作为待处理的数据序列x1(n),x2(n),…,xM(n),n＝0,1,…,N-1，所述的M为阵元通道数，M为大于等于2的正整数，N为检测到的远场窄带信号脉宽长度所对应的采样点个数，取值为2的整数次幂且要求N>2。
第二步：根据给定的阵列参数自动生成等余弦间隔的预引导波束扫描角θk和总的预引导波束个数K，自动生成过程如下，
首先应根据给定的阵列参数计算总的预引导波束个数K，即
Figure BDA0002515294760000021
其中M为均匀线列阵的阵元个数，
Figure BDA0002515294760000022
代表向上取整运算。
然后，依据总的预引导波束个数K自动计算生成预引导波束扫描角θk，即
Figure BDA0002515294760000023
其中arccos[]代表取反余弦运算，k为预引导波束号。
第三步：对已知频率f0的远场窄带信号进行频域波束形成，并求得该已知频率远场窄带信号的各预引导波束扫描角θk所对应的波束功率谱P(θk)，计算过程如下：
首先，计算待处理的均匀线列阵的M通道的远场窄带信号采样数据序列x1(n)，x2(n),…,xM(n),n＝0,1,…,N-1的离散傅里叶变换，其中第m个通道的远场窄带信号采样数据序列xm(n)的离散傅里叶变换Xm(l)的计算方法为：
Figure BDA0002515294760000031
其中m为阵元通道号，m＝1,2,…,M，l为离散频率索引，j表示虚数单位，即
Figure BDA0002515294760000032
然后，依据各预引导波束扫描角θk对待处理的均匀线性阵列M通道的远场窄带信号采样数据序列的离散傅里叶变换X1(l),X2(l),…,XM(l),l＝0,1,…,N-1进行移相处理，其中对第m个通道的远场窄带信号采样数据序列的离散傅里叶变换Xm(l)做与各预引导波束扫描角θk相对应的移相处理结果
Figure BDA0002515294760000033
的计算方法为：
Figure BDA0002515294760000034
其中m为阵元通道号，m＝1,2,…,M，l为离散频率索引，k为预引导波束号，k＝1,2,…,K，d为均匀线列阵的间距，f0为已知的信号频率，c为声音的传播速度。
最后，计算已知频率f0所对应的远场窄带信号各预引导波束扫描角θk所对应的波束功率谱P(θk)，计算过程如下：
Figure BDA0002515294760000035
其中||代表取模值运算，M为均匀线列阵通道数，l0为信号已知频率f0所对应的离散频率索引，其值为：
l0＝round(Nf0/fs) 式(6)
其中round()代表四舍五入运算，fs为设定的信号采样频率。
第四步：搜索波束功率谱P(θk)最大值所对应的预引导波束号kp，搜索过程如下：
根据下式搜索波束功率谱P(θk)最大值所对应的预引导波束号kp：
Figure BDA0002515294760000041
其中
Figure BDA0002515294760000042
表示在1≤k≤K范围内搜索波束功率谱P(θk)的最大值所对应的预引导波束号。
第五步：利用二分数值求解法求解信号来波方向相对偏差
Figure BDA0002515294760000043
过程如下：
首先，计算波束功率谱P(θk)主瓣的峰值和左右相邻的预引导波束所对应的波束功率谱值中较大的一个的比值
Figure BDA0002515294760000044
Figure BDA0002515294760000045
然后，进行参数初始化：设置二分数值求解法的搜索上限a0＝0，搜索下限b0＝0.5，初始化迭代次数q＝1，设置二分精度ε，则二分总迭代次数NB可以表示为：
Figure BDA0002515294760000046
这里
Figure BDA0002515294760000047
表示向上取整运算，ε为大于0且小于0.25的数，log2()表示以2为底的代数运算；
最后，利用二分数值求解法对来波方向相对偏差
Figure BDA0002515294760000048
进行迭代求解，迭代求解过程如下：
(1)令第q次迭代得到的来波方向相对偏差值为
Figure BDA0002515294760000051
(2)计算理论波束功率谱的主瓣峰值和左右相邻的预引导波束所对应的波束功率谱值中较大的一个的比值，计算如下
Figure BDA0002515294760000052
其中,d为均匀线列阵的间距，f0为已知的信号频率，c为声音的传播速度。
(3)重置二分求解法的上下限，即：
如果
Figure BDA0002515294760000053
则令a0＝a0，
Figure BDA0002515294760000054
否则令
Figure BDA0002515294760000055
b0＝b0；
(4)判断是否满足q＝NB，如果满足，转入(5)；否则，令q＝q+1，并返回(1)；
(5)计算求得来波方向相对偏差
Figure BDA0002515294760000056
即
Figure BDA0002515294760000057
第六步：利用
Figure BDA0002515294760000058
估计远场窄带信号的来波方向
Figure BDA0002515294760000059
计算过程如下：
如果1＜kp＜K，则
Figure BDA00025152947600000510
如果kp＝1，则
Figure BDA00025152947600000511
如果kp＝K，则
Figure BDA00025152947600000512
其中
Figure BDA00025152947600000513
为预引导波束扫描角等余弦间隔值。
进一步的，本发明方法的第三步中，第m个通道的远场窄带信号采样数据序列xm(n)的离散傅里叶变换Xm(l)是采用快速傅里叶变换得到的。
进一步的，本发明方法的第五步中，通过利用二分数值求解法对来波方向相对偏差
Figure BDA0002515294760000061
进行迭代求解，二分精度ε设置为0.00001。
本发明方法首先获取M通道阵元数据，然后根据给定的阵列参数生成等余弦间隔的预引导波束扫描角和总的预引导波束个数，接着对已知频率的远场窄带信号进行频域波束形成，并求得各预引导波束扫描角所对应的波束功率谱，搜索最大值对应的预引导波束号，然后利用二分数值求解法求解信号来波方向的相对偏差，最后利用该偏差估计远场窄带信号的来波方向。
有益效果：与现有技术相比，本发明的技术方案具有以下有益技术效果：
1.本发明的估计方法采用理论波束模式模型，因此不存在模型失配的问题。
2.本发明采用二分法数值求解信号来波方向相对偏差从而精确地估计远场窄带信号的来波方向，与抛物线插值算法相比，该算法对相对来波方向偏差具有更强的鲁棒性，并且具有更好的估计性能。
3.本发明仅利用角谱峰值附近的三个样本进行来波方向估计，因而具有与抛物线插值算法相当的计算量。
4.本发明可以推广应用于任意几何阵列的来波方向估计。
附图说明
图1所示为本发明的流程图。
图2所示为M通道阵元接收远场窄带信号示意图。
图3所示为实施例1中仿真M通道阵元的预引导波束扫描角所对应的波束功率谱。
图4所示为实施例2中仿真M通道阵元的预引导波束扫描角所对应的波束功率谱。
具体实施方式
下面结合实施例和说明书附图对本发明作进一步的说明。
如图1所示，一种基于指向性模型数值求解远场窄带信号来波方向估计方法，包括以下步骤：
第一步：获取待处理的均匀线阵列M通道阵元远场窄带信号采样数据序列x1(n)，x2(n),…,xM(n),n＝0,1,…,N-1：从所述阵元基元接收N个采样点的实时采集数据作为待处理的数据序列x1(n),x2(n),…,xM(n),n＝0,1,…,N-1；或从存储器中提取从某一时刻开始接收到的N个信号采样点的数据作为待处理的数据序列x1(n),x2(n),…,xM(n),n＝0,1,…,N-1，所述的M为阵元通道数，M为大于等于2的正整数，N为检测到的远场窄带信号脉宽长度所对应的采样点个数，取值为2的整数次幂且要求N>2。
第二步：根据给定的阵列参数自动生成等余弦间隔的预引导波束扫描角θk和总的预引导波束个数K，自动生成过程如下，
首先应根据给定的阵列参数计算总的预引导波束个数K，即
Figure BDA0002515294760000071
其中M为均匀线列阵的阵元个数，
Figure BDA0002515294760000072
代表向上取整运算；
然后，依据总的预引导波束个数K自动计算生成预引导波束扫描角θk，即
Figure BDA0002515294760000073
其中arccos[]代表取反余弦运算，k为预引导波束序号。
第三步：对已知频率f0的远场窄带信号进行频域波束形成，并求得该已知频率远场窄带信号的各预引导波束扫描角θk所对应的波束功率谱P(θk)，计算过程如下：
首先，计算待处理的均匀线列阵的M通道的远场窄带信号采样数据序列x1(n)，x2(n),…,xM(n),n＝0,1,…,N-1的离散傅里叶变换，其中第m个通道的远场窄带信号采样数据序列xm(n)的离散傅里叶变换Xm(l)的计算方法为：
Figure BDA0002515294760000074
其中m为阵元通道号，m＝1,2,…,M，l为离散频率索引，j表示虚数单位，即
Figure BDA0002515294760000075
然后，依据各预引导波束扫描角θk对待处理的均匀线性阵列M通道的远场窄带信号采样数据序列的离散傅里叶变换X1(l),X2(l),…,XM(l),l＝0,1,…,N-1进行移相处理，其中对第m个通道的远场窄带信号采样数据序列的离散傅里叶变换Xm(l)做与各预引导波束扫描角θk相对应的移相处理结果
Figure BDA0002515294760000081
的计算方法为：
Figure BDA0002515294760000082
其中m为阵元通道号，m＝1,2,…,M，l为离散频率索引，k为预引导波束号，k＝1,2,…,K，d为均匀线列阵的间距，f0为已知的信号频率，c为声音的传播速度；
最后，计算已知频率f0所对应的远场窄带信号各预引导波束扫描角θk所对应的波束功率谱P(θk)，计算过程如下：
Figure BDA0002515294760000083
其中||代表取模值运算，M为均匀线列阵通道数，l0为信号已知频率f0所对应的离散频率索引，其值为：
l0＝round(Nf0/fs) 式(6)
其中round()代表四舍五入运算，fs为设定的信号采样频率。
第四步：搜索波束功率谱P(θk)最大值所对应的预引导波束号kp，搜索过程如下：
根据下式搜索波束功率谱P(θk)最大值所对应的预引导波束号kp：
Figure BDA0002515294760000084
其中
Figure BDA0002515294760000085
表示在1≤k≤K范围内搜索波束功率谱P(θk)的最大值所对应的预引导波束号。
第五步：利用二分数值求解法求解信号来波方向相对偏差
Figure BDA0002515294760000086
过程如下：
首先，计算波束功率谱P(θk)主瓣的峰值和左右相邻的预引导波束所对应的波束功率谱值中较大的一个的比值
Figure BDA0002515294760000087
Figure BDA0002515294760000091
然后，进行参数初始化：设置二分数值求解法的搜索上限a0＝0，搜索下限b0＝0.5，初始化迭代次数q＝1，设置二分精度ε，则二分总迭代次数NB可以表示为：
Figure BDA0002515294760000092
这里
Figure BDA0002515294760000093
表示向上取整运算，ε为大于0且小于0.25的数，log2()表示以2为底的代数运算；
最后，利用二分数值求解法对来波方向相对偏差
Figure BDA0002515294760000094
进行迭代求解，迭代求解过程如下：
(1)令第q次迭代得到的来波方向相对偏差值为
Figure BDA0002515294760000095
(2)计算理论波束功率谱的主瓣峰值和左右相邻的预引导波束所对应的波束功率谱值中较大的一个的比值，计算如下
Figure BDA0002515294760000096
其中d为均匀线列阵的间距，f0为已知的信号频率，c为声音的传播速度。
(3)重置二分求解法的上下限，即：
如果
Figure BDA0002515294760000097
则令a0＝a0，
Figure BDA0002515294760000098
否则令
Figure BDA0002515294760000099
b0＝b0；
(4)判断是否满足q＝NB，如果满足，转入(5)；否则，令q＝q+1，并返回(1)；
(5)计算求得来波方向相对偏差
Figure BDA0002515294760000101
即
Figure BDA0002515294760000102
第六步：利用
Figure BDA0002515294760000103
估计远场窄带信号的来波方向
Figure BDA0002515294760000104
计算过程如下：
如果1＜kp＜K，则
Figure BDA0002515294760000105
如果kp＝1，则
Figure BDA0002515294760000106
如果kp＝K，则
Figure BDA0002515294760000107
其中
Figure BDA0002515294760000108
为预引导波束扫描角等余弦间隔值。
本发明方法首先获取M通道阵元数据，然后根据给定的阵列参数生成等余弦间隔的预引导波束扫描角和总的预引导波束个数，接着对已知频率的远场窄带信号进行频域波束形成，并求得各预引导波束扫描角所对应的波束功率谱，搜索最大值对应的预引导波束号，然后利用二分数值求解法求解信号来波方向的相对偏差，最后利用该偏差估计远场窄带信号的来波方向。
本发明的实施例中，阵元接收远场窄带信号示意图如图2所示，仿真阵元接收远场窄带信号模型为：
Figure BDA0002515294760000109
其中A为均匀线性阵列M通道阵元接收的远场窄带信号的幅度，
Figure BDA00025152947600001010
为阵元接收远场窄带信号的初始相位，N为检测到的远场窄带信号脉宽长度所对应的采样点个数，取值为2的整数次幂且要求N>2，M为总的阵元数目，M为大于等于2的正整数，f0为远场窄带信号频率，fs为采样频率，d为相邻阵元之间的距离，c为远场窄带信号在介质中的传播速度，θT为信号来波方向，即待估计值，wm(n)为第m个阵元接收的相互独立的高斯白噪声，均值为0，方差为σ2且大小由信噪比SNR决定：SNR＝10log10[A2/(2σ2)]。
实施例1
仿真信号参数分别设置为：信号幅度A＝1，初始相位
Figure BDA0002515294760000115
远场窄带信号脉宽长度为1s，接收信号长度为1.032s,采样点个数N＝4096，总的阵元数目M＝21，远场窄带信号频率f0＝1250Hz，采样频率fs＝4000Hz，相邻阵元的间距d＝0.6m，远场窄带信号在介质中的传播速度c＝1500m/s(声音在水中的传播速度)，信号来波方向θT＝50.7035°，SNR＝-6dB。
首先根据给定的阵列参数自动生成总的预引导波束个数K＝25，由此我们可以得到预引导波束扫描角θk，即[0 23.5565 33.5573 41.4096 48.1897 54.3147 60.000065.3757 70.5288 75.5225 80.4059 85.2198 90.0000 94.7802 99.5941 104.4775109.4712 114.6243 120.0000 125.6853 131.8103 138.5904 146.4427 156.4435180.0000]°。
然后对已知频率f0的远场窄带信号进行频域波束形成，并求得该已知频率远场窄带信号的各预引导波束扫描角θk所对应的波束功率谱P(θk)，如图3所示。
接着通过搜索波束功率谱P(θk)最大值，得到所对应的预引导波束号为kp＝5，其对应幅值为
Figure BDA0002515294760000111
相邻两分量幅值分别为
Figure BDA0002515294760000112
Figure BDA0002515294760000113
由此根据式(8)可以计算主瓣的峰值和左右相邻的预引导波束所对应的波束功率谱值中较大的一个的比值为
Figure BDA0002515294760000114
然后设置二分数值求解法的搜索上限a0＝0，搜索下限b0＝0.5，初始化迭代次数q＝1，设置对分精度ε＝0.00001，则根据式(9)可以求出二分总迭代次数为
Figure BDA0002515294760000121
接着利用二分数值求解法对来波方向相对偏差
Figure BDA0002515294760000122
进行迭代求解，得到
Figure BDA0002515294760000123
由于
Figure BDA0002515294760000124
根据式(11)得到来波方向
Figure BDA0002515294760000125
为
Figure BDA0002515294760000126
估计来波方向相对误差为
Figure BDA0002515294760000127
实施例2
仿真信号参数分别设置为：信号幅度A＝1，初始相位
Figure BDA0002515294760000128
远场窄带信号脉宽长度为1s，接收信号长度为1.032s,采样点个数N＝4096，总的阵元数目M＝15，远场窄带信号频率f0＝1250Hz，采样频率fs＝4000Hz，相邻阵元的间距d＝0.6m，远场窄带信号在介质中的传播速度c＝1500m/s(声音在水中的传播速度)，信号来波方向θT＝115.0576°，SNR＝-6dB。
首先根据给定的阵列参数自动生成总的预引导波束个数K＝18，由此我们可以得到预引导波束扫描角θk，即[0 28.0725 40.1192 49.6798 58.0343 65.6843 72.895479.8358 86.6277 93.3723 100.1642 107.1046 114.3157 121.9657 130.3202 139.8808151.9275 180.0000]°。
然后对已知频率f0的远场窄带信号进行频域波束形成，并求得该已知频率远场窄带信号的各预引导波束扫描角θk所对应的波束功率谱P(θk)，结果如图4所示。
接着通过搜索波束功率谱P(θk)最大值，得到所对应的预引导波束号为kp＝13，其对应幅值为
Figure BDA0002515294760000129
相邻两分量幅值分别为
Figure BDA00025152947600001210
Figure BDA00025152947600001211
由此根据式(8)可以计算主瓣的峰值和左右相邻的预引导波束所对应的波束功率谱值中较大的一个的比值为
Figure BDA0002515294760000131
然后设置二分数值求解法的搜索上限a0＝0，搜索下限b0＝0.5，初始化迭代次数q＝1，设置对分精度ε＝0.00001，则根据式(9)可以求出二分总迭代次数为
Figure BDA0002515294760000132
接着利用二分数值求解法对来波方向相对偏差
Figure BDA0002515294760000133
进行迭代求解，得到
Figure BDA0002515294760000134
由于
Figure BDA0002515294760000135
根据式(11)则来波方向
Figure BDA0002515294760000136
为
Figure BDA0002515294760000137
估计来波方向相对误差为
Figure BDA0002515294760000138
从实施例1和例2的结果可以看出，本发明估计方法可以获得良好的估计精度，而且计算简单，计算量小，适用于高精度快速估计多阵元远场窄带信号来波方向的场合。
上述实施例仅是本发明的优选实施方式，应当指出：对于本技术领域的普通技术人员来说，在不脱离本发明原理的前提下，还可以做出若干改进和等同替换，这些对本发明权利要求进行改进和等同替换后的技术方案，均落入本发明的保护范围。

Claims (7)

1.一种基于指向性模型数值求解的远场窄带信号来波方向估计方法，其特征在于，该方法包括以下步骤：
第一步：获取待处理的均匀线阵列M通道阵元远场窄带信号采样数据序列x1(n)，x2(n),…,xM(n),n＝0,1,…,N-1；所述M为阵元通道数，且M为大于等于2的正整数，N为检测到的远场窄带信号脉宽长度所对应的采样点个数，取值为2的整数次幂且要求N>2；
第二步：根据给定的阵列参数自动生成等余弦间隔的预引导波束扫描角θk和总的预引导波束个数K；
第三步：对已知频率f0的远场窄带信号进行频域波束形成，并求得该已知频率远场窄带信号的各预引导波束扫描角θk所对应的波束功率谱P(θk)；
第四步：搜索波束功率谱P(θk)最大值所对应的预引导波束号kp；
第五步：利用二分数值求解法求解信号来波方向相对偏差
Figure FDA0002515294750000011
第六步：利用
Figure FDA0002515294750000012
估计远场窄带信号的来波方向
Figure FDA0002515294750000013
2.按照权利要求1所述的基于指向性模型数值求解的远场窄带信号来波方向估计方法，其特征在于，所述的第一步中，均匀线阵列M通道阵元远场窄带信号采样数据序列x1(n)，x2(n),…,xM(n),n＝0,1,…,N-1的获取方法为：从所述阵列基元接收N个采样点的实时采集数据作为待处理的数据序列x1(n)，x2(n),…,xM(n)；或从存储器中提取从某一时刻开始接收到的N个信号采样点的数据作为待处理的数据序列x1(n)，x2(n),…,xM(n),n＝0,1,…,N-1。
3.按照权利要求2所述的一种基于指向性模型数值求解的远场窄带信号来波方向估计方法，其特征在于，所述的第二步中，根据给定的阵列参数自动生成等余弦间隔的预引导波束扫描角θk和总的预引导波束个数K，方法如下：
(3.1)根据给定的阵列参数计算总的预引导波束个数K，即：
Figure FDA0002515294750000014
其中，M为均匀线列阵的阵元个数，
Figure FDA0002515294750000021
代表向上取整运算；
(3.2)依据总的预引导波束个数K自动计算生成预引导波束扫描角θk，即：
Figure FDA0002515294750000022
其中，arccos[]代表取反余弦运算，k为预引导波束号。
4.按照权利要求3所述的基于指向性模型数值求解的远场窄带信号来波方向估计方法，其特征在于，所述的第三步中，各预引导波束扫描角θk所对应的波束功率谱P(θk)计算过程如下：
(4.1)计算待处理的均匀线列阵的M通道的远场窄带信号采样数据序列x1(n)，x2(n),…,xM(n),n＝0,1,…,N-1的离散傅里叶变换，其中，第m个通道的远场窄带信号采样数据序列xm(n)的离散傅里叶变换Xm(l)的计算方法为：
Figure FDA0002515294750000023
其中，m为阵元通道号，m＝1,2,…,M，l为离散频率索引，j表示虚数单位，即
Figure FDA0002515294750000024
(4.2)依据各预引导波束扫描角θk对待处理的均匀线性阵列M通道的远场窄带信号采样数据序列的离散傅里叶变换X1(l),X2(l),…,XM(l)进行移相处理，其中，对第m个通道的远场窄带信号采样数据序列的离散傅里叶变换Xm(l)做与各预引导波束扫描角θk相对应的移相处理结果
Figure FDA0002515294750000025
的计算方法为：
Figure FDA0002515294750000026
其中，m为阵元通道号，m＝1,2,…,M，l为离散频率索引，k为预引导波束号，k＝1,2,…,K，d为均匀线列阵的间距，f0为已知的信号频率，c为声音的传播速度；
(4.3)计算已知频率f0所对应的远场窄带信号各预引导波束扫描角θk所对应的波束功率谱P(θk)，计算过程如下：
Figure FDA0002515294750000031
其中，| |代表取模值运算，M为均匀线列阵通道数，l0为信号已知频率f0所对应的离散频率索引，其值为：
l0＝round(Nf0/fs) 式(6)
其中，round()代表四舍五入运算，fs为设定的信号采样频率。
5.按照权利要求4所述的基于指向性模型数值求解的远场窄带信号来波方向估计方法，其特征在于，所述的第四步中，搜索波束功率谱P(θk)最大值所对应的预引导波束号kp的方法如下：根据下式搜索波束功率谱P(θk)最大值所对应的预引导波束号kp：
Figure FDA0002515294750000032
其中，
Figure FDA0002515294750000033
表示在1≤k≤K范围内搜索波束功率谱P(θk)的最大值所对应的预引导波束号。
6.按照权利要求5所述的基于指向性模型数值求解的远场窄带信号来波方向估计方法，其特征在于，所述的第五步中，利用二分数值求解法求解信号来波方向相对偏差
Figure FDA0002515294750000034
的方法如下：
(6.1)计算波束功率谱P(θk)主瓣的峰值和左右相邻的预引导波束所对应的波束功率谱值中较大的一个的比值
Figure FDA0002515294750000035
Figure FDA0002515294750000036
(6.2)进行参数初始化：设置二分数值求解法的搜索上限a0，搜索下限b0，初始化迭代次数q＝1，设置二分精度ε，则二分总迭代次数NB可以表示为：
Figure FDA0002515294750000041
其中，
Figure FDA0002515294750000042
表示向上取整运算，ε为大于0且小于0.25的数，log2()表示以2为底的代数运算；
(6.3)利用二分数值求解法对来波方向相对偏差
Figure FDA0002515294750000043
进行迭代求解，迭代求解过程如下：
(1)令第q次迭代得到的来波方向相对偏差值为
Figure FDA0002515294750000044
(2)计算理论波束功率谱的主瓣峰值和左右相邻的预引导波束所对应的波束功率谱值中较大的一个的比值，计算如下：
Figure FDA0002515294750000045
其中，d为均匀线列阵的间距，f0为已知的信号频率，c为声音的传播速度；
(3)重置二分求解法的上下限，即：
如果
Figure FDA0002515294750000046
则令a0＝a0，
Figure FDA0002515294750000047
否则，令
Figure FDA0002515294750000048
b0＝b0；
(4)判断是否满足q＝NB，如果满足，转入(5)；否则，令q＝q+1，并返回(1)；
(5)计算求得来波方向相对偏差
Figure FDA0002515294750000049
即
Figure FDA00025152947500000410
7.按照权利要求6所述的基于指向性模型数值求解的远场窄带信号来波方向估计方法，其特征在于，所述的第六步中，利用
Figure FDA00025152947500000411
估计远场窄带信号的来波方向
Figure FDA00025152947500000412
的方法如下：
如果1＜kp＜K，则
Figure FDA00025152947500000413
如果kp＝1，则
Figure FDA0002515294750000051
如果kp＝K，则
Figure FDA0002515294750000052
其中，
Figure FDA0002515294750000053
为预引导波束扫描角等余弦间隔值。
CN202010474190.6A 2020-05-29 2020-05-29 一种基于指向性模型数值求解的远场窄带信号来波方向估计方法 Active CN111722178B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202010474190.6A CN111722178B (zh) 2020-05-29 2020-05-29 一种基于指向性模型数值求解的远场窄带信号来波方向估计方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202010474190.6A CN111722178B (zh) 2020-05-29 2020-05-29 一种基于指向性模型数值求解的远场窄带信号来波方向估计方法

Publications (2)

Publication Number Publication Date
CN111722178A true CN111722178A (zh) 2020-09-29
CN111722178B CN111722178B (zh) 2022-07-08

Family

ID=72565459

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202010474190.6A Active CN111722178B (zh) 2020-05-29 2020-05-29 一种基于指向性模型数值求解的远场窄带信号来波方向估计方法

Country Status (1)

Country Link
CN (1) CN111722178B (zh)

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113325401A (zh) * 2021-07-06 2021-08-31 东南大学 一种基于线谱相位差解模糊的畸变拖曳线列阵信号重构方法
CN114217265A (zh) * 2021-11-23 2022-03-22 郑州大学 基于最小方差无失真响应的信源波达估计方法及系统

Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN106546949A (zh) * 2016-10-28 2017-03-29 东南大学 一种基于频率预估计的双阵元正弦信号来波方向估计方法
CN106972882A (zh) * 2017-03-01 2017-07-21 浙江大学 基于虚拟域空间功率谱估计的互质阵列自适应波束成形方法
CN107561486A (zh) * 2017-07-03 2018-01-09 西北工业大学 一种基于主动时反的浅海目标波达方向估计方法
CN109188344A (zh) * 2018-08-23 2019-01-11 北京邮电大学 脉冲噪声环境下基于互循环相关music算法信源个数与来波方向角估计方法
CN109639333A (zh) * 2018-12-05 2019-04-16 上海神添实业有限公司 一种基于有效重构协方差矩阵的波束形成方法

Patent Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN106546949A (zh) * 2016-10-28 2017-03-29 东南大学 一种基于频率预估计的双阵元正弦信号来波方向估计方法
CN106972882A (zh) * 2017-03-01 2017-07-21 浙江大学 基于虚拟域空间功率谱估计的互质阵列自适应波束成形方法
CN107561486A (zh) * 2017-07-03 2018-01-09 西北工业大学 一种基于主动时反的浅海目标波达方向估计方法
CN109188344A (zh) * 2018-08-23 2019-01-11 北京邮电大学 脉冲噪声环境下基于互循环相关music算法信源个数与来波方向角估计方法
CN109639333A (zh) * 2018-12-05 2019-04-16 上海神添实业有限公司 一种基于有效重构协方差矩阵的波束形成方法

Cited By (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113325401A (zh) * 2021-07-06 2021-08-31 东南大学 一种基于线谱相位差解模糊的畸变拖曳线列阵信号重构方法
CN113325401B (zh) * 2021-07-06 2024-03-19 东南大学 一种基于线谱相位差解模糊的畸变拖曳线列阵信号重构方法
CN114217265A (zh) * 2021-11-23 2022-03-22 郑州大学 基于最小方差无失真响应的信源波达估计方法及系统

Also Published As

Publication number Publication date
CN111722178B (zh) 2022-07-08

Similar Documents

Publication Publication Date Title
CN111123192B (zh) 一种基于圆形阵列和虚拟扩展的二维doa定位方法
CN108710103B (zh) 基于稀疏阵列的强弱多目标超分辨测向与信源数估计方法
CN109459744B (zh) 一种实现多干扰抑制的稳健自适应波束形成方法
CN104502904B (zh) 一种鱼雷自导波束锐化方法
CN110196407B (zh) 一种基于频率预估的单矢量水听器信号来波方向估计方法
CN111722178B (zh) 一种基于指向性模型数值求解的远场窄带信号来波方向估计方法
CN109581275B (zh) 基于非圆信号和三维正交阵的二维水下doa估计方法和装置
CN111175727B (zh) 一种基于条件波数谱密度的宽带信号方位估计的方法
CN109932679B (zh) 一种传感器列系统最大似然角度分辨率估计方法
CN109407048B (zh) 基于非圆信号和夹角可调阵的水下doa估计方法与装置
CN114137476A (zh) 一种基于极化修正的短波测向方法
CN108957389A (zh) 一种实数域多通道信号目标方位估计方法
CN109884580A (zh) 水下一维doa估计方法和装置
CN109541573A (zh) 一种弯曲水听器阵列的阵元位置校准方法
CN109814065B (zh) 基于相位因子加权的波束形成方法
CN117241375A (zh) 复杂环境中ris辅助基站定位方法、装置、及系统
CN114563760B (zh) 一种基于sca阵型的二阶超波束形成方法、设备及介质
Aboutanios et al. Fast iterative interpolated beamforming for high fidelity single snapshot DOA estimation
CN113126029B (zh) 适用于深海可靠声路径环境的多传感器脉冲声源定位方法
CN111431575B (zh) 基于常规波束形成的来波方向稀疏重构方法
CN110824484B (zh) 一种基于恒模算法的阵元位置估计方法
CN103900682B (zh) 一种时频联合式信号起伏对齐检测方法及系统
CN107238813B (zh) 近场信号源波达方向和波达时间确定方法及装置
CN112114287B (zh) 一种方位观测数据的野值实时剔除方法
CN112649787B (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