CN115276746B - 基于交替方向惩罚的频率一致宽带波束形成方法及系统 - Google Patents
基于交替方向惩罚的频率一致宽带波束形成方法及系统 Download PDFInfo
- Publication number
- CN115276746B CN115276746B CN202210812200.1A CN202210812200A CN115276746B CN 115276746 B CN115276746 B CN 115276746B CN 202210812200 A CN202210812200 A CN 202210812200A CN 115276746 B CN115276746 B CN 115276746B
- Authority
- CN
- China
- Prior art keywords
- dual
- auxiliary
- vector
- constraint
- variable
- 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
- 238000000034 method Methods 0.000 title claims abstract description 117
- 239000013598 vector Substances 0.000 claims abstract description 106
- 238000005457 optimization Methods 0.000 claims abstract description 48
- 239000011159 matrix material Substances 0.000 claims abstract description 17
- 230000009977 dual effect Effects 0.000 claims description 80
- 230000003190 augmentative effect Effects 0.000 claims description 8
- 238000013461 design Methods 0.000 abstract description 6
- 230000006870 function Effects 0.000 description 16
- 238000004088 simulation Methods 0.000 description 16
- 238000010586 diagram Methods 0.000 description 8
- 230000008569 process Effects 0.000 description 6
- 230000008901 benefit Effects 0.000 description 5
- 238000004422 calculation algorithm Methods 0.000 description 4
- 230000000694 effects Effects 0.000 description 4
- 238000010276 construction Methods 0.000 description 3
- 230000001066 destructive effect Effects 0.000 description 3
- 238000004458 analytical method Methods 0.000 description 2
- 238000013459 approach Methods 0.000 description 2
- 230000001427 coherent effect Effects 0.000 description 2
- 230000007423 decrease Effects 0.000 description 2
- 238000012986 modification Methods 0.000 description 2
- 230000004048 modification Effects 0.000 description 2
- NAWXUBYGYWOOIX-SFHVURJKSA-N (2s)-2-[[4-[2-(2,4-diaminoquinazolin-6-yl)ethyl]benzoyl]amino]-4-methylidenepentanedioic acid Chemical compound C1=CC2=NC(N)=NC(N)=C2C=C1CCC1=CC=C(C(=O)N[C@@H](CC(=C)C(O)=O)C(O)=O)C=C1 NAWXUBYGYWOOIX-SFHVURJKSA-N 0.000 description 1
- KLXQAXYSOJNJRI-KVTDHHQDSA-N (2s,3s,4r,5r)-5-amino-2,3,4,6-tetrahydroxyhexanal Chemical compound OC[C@@H](N)[C@@H](O)[C@H](O)[C@H](O)C=O KLXQAXYSOJNJRI-KVTDHHQDSA-N 0.000 description 1
- 101000651373 Homo sapiens Serine palmitoyltransferase small subunit B Proteins 0.000 description 1
- 102100027676 Serine palmitoyltransferase small subunit B Human genes 0.000 description 1
- 239000000654 additive Substances 0.000 description 1
- 230000000996 additive effect Effects 0.000 description 1
- 230000003321 amplification Effects 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 238000004364 calculation method Methods 0.000 description 1
- 230000008859 change Effects 0.000 description 1
- 230000003247 decreasing effect Effects 0.000 description 1
- 230000001419 dependent effect Effects 0.000 description 1
- 239000006185 dispersion Substances 0.000 description 1
- 238000012804 iterative process Methods 0.000 description 1
- 238000003199 nucleic acid amplification method Methods 0.000 description 1
- 238000012545 processing Methods 0.000 description 1
- 238000003672 processing method Methods 0.000 description 1
- 238000011160 research Methods 0.000 description 1
- 230000004044 response Effects 0.000 description 1
- 238000005070 sampling Methods 0.000 description 1
- 230000035945 sensitivity Effects 0.000 description 1
- 230000002123 temporal effect Effects 0.000 description 1
- 238000012360 testing method Methods 0.000 description 1
- 230000000007 visual effect Effects 0.000 description 1
Images
Classifications
-
- H—ELECTRICITY
- H04—ELECTRIC COMMUNICATION TECHNIQUE
- H04B—TRANSMISSION
- H04B7/00—Radio transmission systems, i.e. using radiation field
- H04B7/02—Diversity systems; Multi-antenna system, i.e. transmission or reception using multiple antennas
- H04B7/04—Diversity systems; Multi-antenna system, i.e. transmission or reception using multiple antennas using two or more spaced independent antennas
- H04B7/06—Diversity systems; Multi-antenna system, i.e. transmission or reception using multiple antennas using two or more spaced independent antennas at the transmitting station
- H04B7/0613—Diversity systems; Multi-antenna system, i.e. transmission or reception using multiple antennas using two or more spaced independent antennas at the transmitting station using simultaneous transmission
- H04B7/0615—Diversity systems; Multi-antenna system, i.e. transmission or reception using multiple antennas using two or more spaced independent antennas at the transmitting station using simultaneous transmission of weighted versions of same signal
- H04B7/0617—Diversity systems; Multi-antenna system, i.e. transmission or reception using multiple antennas using two or more spaced independent antennas at the transmitting station using simultaneous transmission of weighted versions of same signal for beam forming
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
- G06F17/10—Complex mathematical operations
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
- G06F17/10—Complex mathematical operations
- G06F17/15—Correlation function computation including computation of convolution operations
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/20—Design optimisation, verification or simulation
-
- 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
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- General Physics & Mathematics (AREA)
- Data Mining & Analysis (AREA)
- Mathematical Physics (AREA)
- General Engineering & Computer Science (AREA)
- Computational Mathematics (AREA)
- Mathematical Analysis (AREA)
- Mathematical Optimization (AREA)
- Pure & Applied Mathematics (AREA)
- Databases & Information Systems (AREA)
- Software Systems (AREA)
- Algebra (AREA)
- Computer Hardware Design (AREA)
- Evolutionary Computation (AREA)
- Geometry (AREA)
- Computing Systems (AREA)
- Computer Networks & Wireless Communication (AREA)
- Signal Processing (AREA)
- Variable-Direction Aerials And Aerial Arrays (AREA)
- Radar Systems Or Details Thereof (AREA)
Abstract
本发明公开了基于交替方向惩罚的频率一致宽带波束形成方法及系统,该方法包括:步骤一,以hH(ωn)Γ(ωn)h(ωn)最小为目标函数,结合无失真约束、主瓣约束、副瓣约束以及白噪声增益约束构建优化问题;其中,h(ωn)表示角频率ωn对应的权值向量,Γ(ωn)表示伪相干矩阵;步骤二,利用交替方向惩罚法求解优化问题获得最优权值向量;步骤三,将最优的权值向量与导向矢量相乘结果作为频率一致宽带波束形成。本发明设计了频率一致波束形成器,可以确保近似频率一致以及频率相关的白噪声增益值;且对交替方向法进行改进,可以更高效的求解优化问题。
Description
技术领域
本发明属于波束形成技术领域,尤其涉及基于交替方向惩罚的频率一致宽带波束形成方法及系统。
背景技术
宽带波束形成是声学领域的重要研究内容,包括频率不变波束形成(FIB)。频率不变波束形成技术可在接受目标信号范围内形成一定宽度的波束,以抑制信号失真,在实际中可用于水声中波场的采样。
为满足实际需要,超指向波束形成技术被开发出来。当线性阵列的阵元间距趋于0时,其方向性因子(DF)在端射方向达到最大值,此时紧凑的孔径易于嵌入到其他器件中,同时提供高的方向增益。然而,超指向波束形成器在低频时仍然对白噪声和阵列不匹配很敏感。为解决该问题,许多设计波束形成器的优化问题被开发出来,优化问题常用的求解方法之一是内点法,但内点法往往会耗费大量的计算时间。
近年来,交替方向乘子法(ADMM)在求解优化问题方面也取得了一定的成就,可以实现高效求解非凸非线性优化问题,但存在不收敛问题,且运算过程也有待简化。
发明内容
本发明的目的是提供基于交替方向惩罚的频率一致宽带波束形成方法及系统。
本发明的技术方案如下:
基于交替方向惩罚的频率一致宽带波束形成方法,包括:
步骤一,以hH(ωn)Γ(ωn)h(ωn)最小为目标函数,结合无失真约束、主瓣约束、副瓣约束以及白噪声增益约束构建优化问题;其中,h(ωn)表示角频率ωn对应的权值向量,Γ(ωn)表示伪相干矩阵;
步骤二,利用交替方向惩罚法求解优化问题获得最优权值向量,该步骤进一步包括:
2.1引入似实值助变量和对偶辅助向量重新表达优化问题;
2.2引入各线性约束的缩放对偶变量构造优化问题的增广拉格朗日函数;
2.3对权值向量、似实值助变量、对偶辅助向量、惩罚因子、缩放对偶变量进行交叉迭代求解;
2.4重复子步骤2.3直至达到收敛条件,此时的权值向量即最优的权值向量;
步骤三,将最优的权值向量与导向矢量相乘结果作为频率一致宽带波束形成。
进一步的,步骤一所构建的优化问题如下:
其中:
h(ωn)表示角频率ωn对应的权值向量,Γ(ωn)表示伪相干矩阵;
hH(ωn)d(ωn,cosθd)=1为无失真约束,d(ωn,cosθd)为角频率ωn对应的导向矢量,θd为期望信号的入射方向角;
hH(ωn)h(ωn)≤γ(ωn)为白噪声增益约束,γ(ωn)为角频率ωn处白噪声增益的上限;
以及子步骤2.1中重新表达的优化问题如下:
其中:s1(ωn)和s2(ωn)为引入的似实值助变量,z(ωn)为引入的对偶辅助向量。
进一步的,子步骤2.2中构造的增广拉格朗日函数如下:
其中:ρ1、ρ2、ρ3、ρ4均为自定义的惩罚因子,也即迭代步长,ρ1>0,ρ2>0,ρ3>0,ρ4>0;μ(ωn)为拉格朗日乘子;
μ(ωn)为线性约束hH(ωn)d(ωn,cosθd)=1对应的缩放对偶变量;
v(ωn)为线性约束z(ωn)=h(ωn)对应的缩放对偶变量;
基于上述构造的增广拉格朗日函数,子步骤2.3进一步包括:
2.3a初始化权值向量、似实值助变量、对偶辅助向量、惩罚因子、缩放对偶变量;
2.3b基于似实值助变量、对偶辅助向量、惩罚因子、缩放对偶变量的当前值,通过求解无约束优化问题,对权值向量进行更新;
2.3c基于权值向量、似实值助变量、惩罚因子、缩放对偶变量的当前值,仅在与对偶辅助向量相关的约束条件下求解优化问题,对偶辅助向量进行更新;
2.3d基于权值向量、对偶辅助向量、惩罚因子、缩放对偶变量的当前值,仅在与似实值助变量相关的约束条件下求解优化问题,对似实值助变量进行更新;
2.3e基于权值向量、似实值助变量、对偶辅助向量、惩罚因子的当前值,更新缩放对偶变量;
2.4f基于权值向量、、似实值助变量、对偶辅助向量、缩放对偶变量的当前值,更新惩罚因子。
其中,h(k+1)(ωn)为更新后的权值向量;Ξ1和Ξ2如下:
更新后的对偶辅助变量zk+1(ωn)如下:
进一步的,子步骤2.3d中,先基于权值向量、第二似实值助变量、对偶辅助向量、惩罚因子、缩放对偶变量的当前值,仅在与第一似实值助变量相关的约束条件下求解优化问题,对第一似实值助变量进行更新;
再基于权值向量、第一似实值助变量、对偶辅助向量、惩罚因子、缩放对偶变量的当前值,仅在与第二似实值助变量相关的约束条件下求解优化问题,对第二似实值助变量进行更新。
进一步的,子步骤2.4中收敛条件为:
本发明提供的基于交替方向惩罚的频率一致宽带波束形成系统,包括:
第一模块,用来以hH(ωn)Γ(ωn)h(ωn)最小为目标函数,结合无失真约束、主瓣约束、副瓣约束以及白噪声增益约束构建优化问题;其中,h(ωn)表示角频率ωn对应的权值向量,Γ(ωn)表示伪相干矩阵;
第二模块,用来利用交替方向惩罚法求解优化问题获得最优权值向量;
所述第二模块进一步包括子模块:
第一子模块,用来引入似实值助变量和对偶辅助向量重新表达优化问题;
第二子模块,用来引入各线性约束的缩放对偶变量构造优化问题的增广拉格朗日函数;
第三子模块,用来对权值向量、似实值助变量、对偶辅助向量、惩罚因子、缩放对偶变量进行交叉迭代求解;
第四子模块,用来判断重复第三子模块直至达到收敛条件,此时的权值向量即最优的权值向量;
第三模块,用来将最优的权值向量与导向矢量相乘结果作为频率一致宽带波束形成。
和现有技术相比,本发明具有如下优点和有益效果:
本发明设计了频率一致波束形成器,可以确保近似频率一致以及频率相关的白噪声增益值。目前优化问题求解方法普遍存在效率低且非收敛的问题,加之在不同频率上搜索惩罚参数也存在一定难度,本发明对交替方向法进行改进,提出了一种交替方向惩罚法来求解优化问题。该方法利用残差和惩罚因子间的关系来减少迭代次数,可以更高效的求解优化问题,尤其是在求解具大量阵元优化问题时,求解效率更是得到极大提高。
附图说明
图1为本发明方法流程示意图;
图8为本发明方法在M=8,δ=1cm仿真条件下的宽带波束图和在f=2500Hz下的极坐标波束图,其中,图(a)为宽带波束度,图(b)为极坐标波束图;
图9为非相干法在M=8,δ=1cm仿真条件下的宽带波束图和在f=2500Hz下的极坐标波束图,其中,图(a)为宽带波束度,图(b)为极坐标波束图;
图10为MN法在M=8,δ=1cm仿真条件下的宽带波束图和在f=2500Hz下的极坐标波束图,其中,图(a)为宽带波束度,图(b)为极坐标波束图;
图11为ICE法在M=8,δ=1cm仿真条件下的宽带波束图和在f=2500Hz下的极坐标波束图,其中,图(a)为宽带波束度,图(b)为极坐标波束图;
图12为本发明方法、非相干法、MN法、ICE法的DF曲线;
图13为本发明方法、非相干法、MN法、ICE法的WNG曲线;
图14为本发明方法在M=8,δ=1cm仿真条件下α为0.2时的宽带波束图和在f=2500Hz下的极坐标波束图,其中,图(a)为宽带波束度,图(b)为极坐标波束图;
图15为本发明方法在M=8,δ=1cm仿真条件下α为0.45时的宽带波束图和在f=2500Hz下的极坐标波束图,其中,图(a)为宽带波束度,图(b)为极坐标波束图;
图16为本发明方法在M=8,δ=1cm仿真条件下α为0.6时的宽带波束图和在f=2500Hz下的极坐标波束图,其中,图(a)为宽带波束度,图(b)为极坐标波束图;
图17为本发明方法在M=8,δ=1cm仿真条件下α为0.8时的宽带波束图和在f=2500Hz下的极坐标波束图,其中,图(a)为宽带波束度,图(b)为极坐标波束图;
图18为α取0.2、0.45、0.6、0.8时本发明方法的DF;
图19为α取0.2、0.45、0.6、0.8时本发明方法的WNG曲线。
具体实施方式
为了使本发明的目的、技术方案及优点更加清楚明白,以下结合附图对本发明具体实施方式进行详细说明。应当理解,此处所描述的具体实施方式仅仅用以解释本发明,并不限定本发明。此外,下面所描述的具体实施方式中所涉及到的技术特征只要彼此之间未构成冲突就可以相互组合。
本文中所涉及字符参数含义如下:
小写的粗体符号表示向量,大写的粗体符号表示矩阵,(大写)。为虚数单位;(·)T、(·)H、(·)*分别表示转置算子、共轭转置算子、复共轭算子;vec(·)表示向量化运算符;IM表示大小M×M的单位矩阵;||·||2表示二范数;|·|表示模值;O(·)表示算法复杂度的朗道符号(Landau符号)。
为便于理解,下面将对本发明所基于的技术原理及所采用技术方案分别详细说明。
一、信号模型构建
考虑基于均匀间隔的线性阵列的宽带波束形成结构,远场源信号(平面波)以声速c=340m/s入射到由全方位传感器组成的阵列上。两个连续传感器间的距离等于δ,远场源信号到阵列的入射方位角采用θ表示,导向矢量d(ω,cosθ)由式(1)给出:
式(1)中:
ω表示角频率,ω=2πf,f表示时间频率,f>0,则波长λ=c/f;
τ0=δ/c是入射方位角θ=0°处两个连续传感器间的时间差;
e表示自然对数,M表示阵列中传感器数量,即阵列中阵元数量。
式(2)中;
假设期望信号从入射方位角θd传播,则观测信号向量y(ω)为:
式(3)中:Ym(ω)表示第m个传感器的观测信号;d(ω,cosθd)表示期望信号源所在方向的导向矢量;X(ω)表示期望信号;vN(ω)表示加性噪声信号矢量。
通常波束形成是通过对每个传感器信号应用时间滤波器,并对滤波后的信号求和来完成。在频域中,相当于向每个传感器的输出添加一个复数权重,并对整个阵元直径求和,即:
式(4)中,Z(ω)表示波束形成器输出信号;Xfd(ω)表示经过滤波器的所需信号,Vrn(ω)表示残余噪声。
二、优化问题构建
式(5)中,Ω和Θ分别表示感兴趣的角频率范围和入射角度范围;εLS为正公差参数,表示合成波束图和参考波束图间的误差。
虽然式(5)中合成波束图考虑到所有的入射角度和角频率,但不能考虑不同角频率下获得频率一致波束的能力差异,主瓣区域的频率不变性比副瓣区域的频率不变性更重要,这意味着必须优先考虑以确保主瓣的频率不变性。因此,统一离散角频率和入射角度,引入I个方向和J个角频率单元,将主瓣和副瓣分别划分为L和S个角度。式(5)重写为:
式(6)中,上标m表示主瓣,上标s表示副瓣,Θm表示主瓣区域,Θs表示副瓣区域,C是表示微分dωdθ的正常数;
Dm(ωn)表示大小M×L的导向矩阵,即角频率ωn下覆盖主瓣区域的导向矢量所组成的导向矩阵,其表示见式(7):
Dm(ωn)=[d(ωn,cosθ1),d(ωn,cosθ2),...,d(ωn,cosθL)] (7)
式(8)中,ωn表示第n个频率,h(ωn)则为角频率ωn对应的复数权值向量,d(ωn,cosθd)为角频率ωn对应的导向矢量。
无失真约束意味:任何到达的信号将沿着导向矢量d(ωj,cosθd)在不失真情况下通过波束形成器。
式(9)~(10)中,ε1(ωn)和ε2(ωn)分别为主瓣和副瓣对应的正误差参数。
由于频率不变性决定了各种性能指标,且主瓣比副瓣的影响更大,根据式(5)和(6),建立以下关系:
它保证了阵列的鲁棒性,γ(ωn)表示在角频率ωn处白噪声增益的最大值上限。
对于上述优化问题,一般处理方法采用内点法,但优化问题中所有约束相互耦合,因为都为向量h(ωn)的函数,采用内点法的计算复杂度为O(M3.5(L+S)),这会耗费大量的计算时间,因此需要采用一种快速收敛方法来解决该问题。
三、基于交替方向惩罚的频率一致波束形成方法
利用目标函数和线性约束以缩放形式来构造优化问题(16)的增广拉格朗日函数:
式(17)中:ρ1、ρ2、ρ3、ρ4均为自定义的惩罚因子,也即迭代步长,ρ1>0,ρ2>0,ρ3>0,ρ4>0;μ(ωn)为拉格朗日乘子;
μ(ωn)为线性约束hH(ωn)d(ωn,cosθd)=1对应的缩放对偶变量;
v(ωn)为线性约束z(ωn)=h(ωn)对应的缩放对偶变量;
更新在每次迭代中进行,而不是在双重上升法中直到内部循环收敛的多次更新。更重要的是,本发明交替方向惩罚在ADMM框架上完成,不仅限于以上常见的优点。交替方向惩罚法基于原始残差更新惩罚因子和缩放对偶变量,进一步提高收敛能力,这种省时的实现策略可以确保交替方向惩罚法有更高概率找到比ADMM和内点法更优的可行点。
本发明利用交替方向惩罚法(ADMP)求解优化问题的过程如下:
S200:通过求解式(18)所示的无约束优化问题,来更新h(ωn):
每次迭代中都有h(ωn)闭合解,h(ωn)更新公式如下:
其中:
将式(17)代入式(23)并移除无关项,式(23)改写为:
采用式(24)更新对偶辅助变量z(ωn):
步骤400:使用hk+1(ωn),zk+1(ωn),s2 k(ωn),μk(ωn),vk(ωn),wk(ωn),rk(ωn)解决式(25)所示优化问题并更新s1(ωn):
将式(17)插入(25)并移除无关项,则优化问题被改写为:
s1(ωn)的更新公式如下:
同理可得s2(ωn)的更新公式如下:
步骤500:采用式(29)~(30)更新缩放对偶变量μ(ωn)、v(ωn)、w(ωn)、r(ωn):
式(29)~(32)中涉及变量定义如下:
其中,η2、η3、η4为足够大的正数。
步骤600:采用式(34)更新惩罚因子:
参数δa1、δa2的值均应接近1,该值与收敛速度及相应的残值有关,而惩罚因子的更新与原始残差有关。具体而言,如果不随迭代次数减少,则使用更大的惩罚参数/>使/>接近零,从而以高概率找到行点;否则,/>保持不变。
上标k和k+1分别表示该变量为迭代更新前后的值。
分析本发明所采用算法的计算复杂度,每次迭代中对h(ωn)、Ξ1、Ξ2的更新,计算复杂度分别为O(M2L)、O(M3)和O(ML),步骤1的计算复杂度共计O(ML+M2L+M3);步骤2的计算仅涉及简单的向量运算,其复杂度为O(M);步骤3是传统ADMM算法的一部分,与步骤1~2相同,每次迭代的复杂度为O(ML+MS);步骤4与步骤5相似,每次迭代的复杂度为O(2M+ML+MS)。因此,本发明算法每次迭代的计算复杂度共计O(M3+M2L+3ML+4MS+5M)。
四、仿真实验
本仿真实验在均匀线阵比较本发明方法、非相干法、最小范数法(简记为MN法)、不相关拓展法(简记为ICE法)的性能,非相干法是以hH(ωn)h(ωn)最小为目标函数,以为约束条件来构建优化问题,并采用内点法求解。MN法是基于N阶差分麦克风的最小范数设计,其将(差分麦克风)DMA的阶数与传感器数量不关联化,并在固定阶数下增加阵元数。ICE法与非相干法类似,唯一区别在于ICE法中误差参数ε1(ωn)、ε2(ωn)和阵元间距δ相对增加。
将这四种方法用于在flow=200Hz和fhigh=4480Hz间的频率范围内,设计频率一致的宽带波束图。假设用于相应时域接收信号的窗口分析的持续时间T=25毫秒,频率分辨率f=1/T=40Hz。因此,频点数量最后,考虑为均匀线阵,使用8个全向麦克风,阵元间距/>/>
式(25)中,aN,n是确定N阶差分麦克风不同方向图的实归一化系数值,n=0,1,...,N;aN是长度向量,aN=[aN,0,aN,1,...aN,N]T。
该三阶锐心形波束图可表示为B3,Hd(θ),B3,Hd(θ)为期望的波束图:
B3,Hd(θ)=-0.14-0.57cosθ+0.57cos2θ+1.15cos3θ (26)
建立矩阵Dm(ωn)和Ds(ωn),令离散角Δθ=1°,ρL=θL=60°,ρs=θs=120°,即主瓣区域是[-60°,-60°],旁瓣区域是[-180°,-61°]∪[61°,180°]∪[61°,180°]。
对于本发明方法,设置参数为:ε1(ωn)=0.1,ε2(ωn)=0.2,α(ωn)=0.4,Δr1=Δr2=Δr3=Δr4=10-6,Δ1=10-6,Δ2=10-3,δ11=0.9985,δ21=0.9995,δ31=0.9995,δ41=0.9995,δ12=1.025,δ22=1.005,δ32=1.005,δ42=1.005。对于非相干法,Δra、Δ1、Δ2取值同本发明,但δ11=0.95,δ12=1.05。
求解时目标函数在一定范围收敛。为测试迭代次数对收敛性的影响,令采用本发明方法进行求解,仿真条件为:f=200Hz,M=8,δ=1cm。不同收敛范围下的收敛结果见图2,原始残差Er1(ωn)与迭代次数的关系见图3,原始残差Er3(ωn)与迭代次数的关系见图4,原始残差Er4(ωn)与迭代次数的关系见图5,对偶残差Eh(ωn)与迭代次数的关系见图6,仿真波束图见图7。
图3~6给出了残差随迭代次数的变化规律,从图中可以看出,残差岁迭代次数而降低,且四个残差最终都趋于0,得到了收敛。更准确来说,在一定迭代次数后,原始残差变为0,而对偶残差接近0。通过比较对收敛速度的影响,从图2~6可以看出,少量改变/>对收敛速度有显著影响,而/>值的改变对收敛速度只有少量影响,在迭代过程中/>比/>对目标函数收敛速度影响更大,于是令/>但/>
图2给出了七种下的目标值DF,从图2可以看出,随着迭代次数增加,均可收敛到11.7dB~11.9dB。因为本发明方法所构建的优化问题是在整个宽带上,处理问题的效率很重要,因此仅考虑目标函数的值是不够的,更希望能在目标函数值和收敛速度上折中。目标函数值与惩罚因子/>没有直观的关系,这也是由于惩罚系数也随着迭代次数变化,这点不同于传统的ADMM。见图7,图中虚线图形代表本发明方法所求波束图,实线图像代表期望波束图;从图7可以观察到,在特定的惩罚参数/>下设计的波束图,低频时非常接近期望波束图,这意味着本发明可以很好地求解所构建的优化问题。
见图8~11,分别为上述四种方法在M=8,δ=1cm仿真条件下的宽带波束图和在f=2500Hz下的极坐标波束图。其中图8(b)、图9(b)、图10(b)、图11(b)展示了特定频率f=2500Hz下极坐标波束图(见虚线图形)与期望波束图(见实线图形)(即,三阶锐心形波束图)的比较。本发明方法、非相干法、MN法的阵元线阵布局相同,ICE法的误差参数和阵元间距在非相关方法的基础上发生变化。从图中可以看出,MN法所得波束图与期望波束图最相似,而ICE法所得波束图具有更高旁瓣,本发明利用ICE法可找到一个有效的白噪声输出功率上限,但会牺牲波束图频率一致,以便能在特定范围内获得效果较好的频率一致波束图。对于有着更好白噪声增益性能的非相关法,其所得波束图的波瓣与期望波束图波瓣几乎相似,但代价是旁瓣和主瓣间的误差相对较大。
见图8(a)、图9(a)、图10(a)、图11(a),本仿真实验中,期望波束图在200Hz和2500Hz间的频率范围内具有三个零点,具体而言,在200Hz和1000Hz间,波束图与期望波束图案几乎相似,第二个和第三个零点的角度逐渐增大,第三个零点从154°增加到180°直至消失。总体而言,随着频率增加,零点的数量大约从f=2500Hz开始减少到两个,从深到浅变化。对于本发明方法和非相干法,所得宽带波束图几乎频率一致,且在主瓣区域误差ε1(ωn)较小,在旁瓣区域误差ε2(ωn)较大。MN法实现了最完美的频率一致波束图,波束方向图在所有频率的主瓣和旁瓣范围内都是一致的;对于ICE法,对于大部分频率都有两个零点,第二个零点在3000Hz和4480Hz不明显,频率一致性较差。
图12~13为四种方法在M=8,δ=1cm仿真条件下的DF和WNG曲线,通过与期望波束图的DF和WNG曲线比较,本发明方法和不相干法具有与期望波束图最相似的结果。
本发明方法具有比不相干法更高效的运算速度,同时在64GB内存的AMD Ryzen 93900X 12Core上运行本发明方法和不相干法,在不包括参数搜索过程的情况下,本发明方法大约需要15秒运行时间,而不相干法则运行了将近90秒,几乎是本发明方法的六倍。虽然MN法可获得频率一致的相对最优DF,这意味着它的波束图频率一致性完美,但会遭受低频白噪声的放大,对于较少的线性阵列阵元,DF几乎恒定,最高可达8000Hz;但随着M的增加,DF从大约f=2000Hz开始不断减小;但在高频下,本发明方法可获得更高的WNG。
图14~图17展示了在折中参数α为0.2、0.45、0.6、0.8时本发明方法所获得的宽带波束图和f=2500Hz下的极坐标波束图,在极坐标波束图中,实现图形为期望波束图。图18~图19展示了α为0.2、0.45、0.6、0.8时本发明方法的DF和WNG曲线。从图14~19可以看出,本发明FI超指向性波束形成器在WNG和DF间取得了很好的折中,特别是α为0.45时。随着α增加,以牺牲WNG为代价,DF得到改善,并且波束图零点更清晰。
本发明所提出的频率不变波束形成器,所形成波束图主瓣区域几乎频率不变,同时还可以实现所需的WNG值,并具有一定的灵活性以最大化DF。与不相关法相比,本发明在显著降低计算复杂度情况下还实现了良好的性能,并且优化方法可用于所有方向和感兴趣频率。
上述实施例所述是用以具体说明本发明,文中虽通过特定的术语进行说明,但不能以此限定本发明的保护范围,熟悉此技术领域的人士可在了解本发明的精神与原则后对其进行变更或修改而达到等效目的,而此等效变更和修改,皆应涵盖于权利要求范围所界定范畴内。
Claims (5)
1.基于交替方向惩罚的频率一致宽带波束形成方法,其特征是,包括:
步骤一,以hH(ωn)Γ(ωn)h(ωn)最小为目标函数,结合无失真约束、主瓣约束、副瓣约束以及白噪声增益约束构建优化问题;其中,h(ωn)表示角频率ωn对应的权值向量,Γ(ωn)表示伪相干矩阵;
步骤二,利用交替方向惩罚法求解优化问题获得最优权值向量,该步骤进一步包括:
2.1引入似实值助变量和对偶辅助向量重新表达优化问题;
2.2引入各线性约束的缩放对偶变量构造优化问题的增广拉格朗日函数;
2.3对权值向量、似实值助变量、对偶辅助向量、惩罚因子、缩放对偶变量进行交叉迭代求解;
2.4重复子步骤2.3直至达到收敛条件,此时的权值向量即最优的权值向量;
步骤三,将最优的权值向量与导向矢量相乘结果作为频率一致宽带波束形成;
步骤一所构建的优化问题如下:
其中:
h(ωn)表示角频率ωn对应的权值向量,Γ(ωn)表示伪相干矩阵;
hH(ωn)d(ωn,cosθd)=1为无失真约束,d(ωn,cosθd)为角频率ωn对应的导向矢量,θd为期望信号的入射方向角;
hH(ωn)h(ωn)≤γ(ωn)为白噪声增益约束,γ(ωn)为角频率ωn处白噪声增益的上限;
以及子步骤2.1中重新表达的优化问题如下:
其中:s1(ωn)和s2(ωn)为引入的似实值助变量,z(ωn)为引入的对偶辅助向量;
子步骤2.2中构造的增广拉格朗日函数如下:
其中:ρ1、ρ2、ρ3、ρ4均为自定义的惩罚因子,也即迭代步长,ρ1>0,ρ2>0,ρ3>0,ρ4>0;μ(ωn)为拉格朗日乘子;
μ(ωn)为线性约束hH(ωn)d(ωn,cosθd)=1对应的缩放对偶变量;
v(ωn)为线性约束z(ωn)=h(ωn)对应的缩放对偶变量;
基于上述构造的增广拉格朗日函数,子步骤2.3进一步包括:
2.3a初始化权值向量、似实值助变量、对偶辅助向量、惩罚因子、缩放对偶变量;
2.3b基于似实值助变量、对偶辅助向量、惩罚因子、缩放对偶变量的当前值,通过求解无约束优化问题,对权值向量进行更新;
2.3c基于权值向量、似实值助变量、惩罚因子、缩放对偶变量的当前值,仅在与对偶辅助向量相关的约束条件下求解优化问题,对偶辅助向量进行更新;
2.3d基于权值向量、对偶辅助向量、惩罚因子、缩放对偶变量的当前值,仅在与似实值助变量相关的约束条件下求解优化问题,对似实值助变量进行更新;
2.3e基于权值向量、似实值助变量、对偶辅助向量、惩罚因子的当前值,更新缩放对偶变量;
2.3f基于权值向量、似实值助变量、对偶辅助向量、缩放对偶变量的当前值,更新惩罚因子;
子步骤2.4中收敛条件为:
2.如权利要求1所述的基于交替方向惩罚的频率一致宽带波束形成方法,其特征是:
其中,h(k+1)(ωn)为更新后的权值向量;Ξ1和Ξ2如下:
4.如权利要求1所述的基于交替方向惩罚的频率一致宽带波束形成方法,其特征是:
子步骤2.3d中,先基于权值向量、第二似实值助变量、对偶辅助向量、惩罚因子、缩放对偶变量的当前值,仅在与第一似实值助变量相关的约束条件下求解优化问题,对第一似实值助变量进行更新;
再基于权值向量、第一似实值助变量、对偶辅助向量、惩罚因子、缩放对偶变量的当前值,仅在与第二似实值助变量相关的约束条件下求解优化问题,对第二似实值助变量进行更新。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202210812200.1A CN115276746B (zh) | 2022-07-12 | 2022-07-12 | 基于交替方向惩罚的频率一致宽带波束形成方法及系统 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202210812200.1A CN115276746B (zh) | 2022-07-12 | 2022-07-12 | 基于交替方向惩罚的频率一致宽带波束形成方法及系统 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN115276746A CN115276746A (zh) | 2022-11-01 |
CN115276746B true CN115276746B (zh) | 2023-05-30 |
Family
ID=83765465
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN202210812200.1A Active CN115276746B (zh) | 2022-07-12 | 2022-07-12 | 基于交替方向惩罚的频率一致宽带波束形成方法及系统 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN115276746B (zh) |
Families Citing this family (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN117406173B (zh) * | 2023-12-11 | 2024-03-15 | 湘江实验室 | 一种声源检测用波束形成方法及装置 |
Citations (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN110932734A (zh) * | 2019-11-14 | 2020-03-27 | 浙江大学 | 一种基于交替方向乘子法的深度学习信道译码方法 |
Family Cites Families (7)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US9753892B2 (en) * | 2014-02-20 | 2017-09-05 | Mitsubishi Electric Research Laboratories, Inc. | Method for solving quadratic programs for convex sets with linear equalities by an alternating direction method of multipliers with optimized step sizes |
CN106682405B (zh) * | 2016-12-14 | 2019-10-18 | 西北工业大学 | 基于凸优化的低旁瓣波束图综合设计方法 |
CN108735228B (zh) * | 2017-04-20 | 2023-11-07 | 斯达克实验室公司 | 语音波束形成方法及系统 |
CN108809397B (zh) * | 2018-06-27 | 2020-06-30 | 东南大学 | 多天线系统中高功效数模混合波束成形方法、装置及设备 |
CN111551923B (zh) * | 2020-05-27 | 2022-11-04 | 电子科技大学 | 一种多约束下的均匀线性阵列低旁瓣波束形成优化方法 |
CN113704998B (zh) * | 2021-08-27 | 2023-03-03 | 电子科技大学 | 一种基于交替方向乘子法的波束形成方法 |
CN115130331B (zh) * | 2022-08-30 | 2022-12-13 | 湖北工业大学 | 基于稀疏阵列的鲁棒频率不变波束形成方法 |
-
2022
- 2022-07-12 CN CN202210812200.1A patent/CN115276746B/zh active Active
Patent Citations (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN110932734A (zh) * | 2019-11-14 | 2020-03-27 | 浙江大学 | 一种基于交替方向乘子法的深度学习信道译码方法 |
Also Published As
Publication number | Publication date |
---|---|
CN115276746A (zh) | 2022-11-01 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN106788653B (zh) | 一种基于协方差矩阵重构的自适应波束形成方法 | |
EP1429581B1 (en) | Method of broadband constant directivity beamforming for non linear and non axi-symmetric sensor arrays embedded in a obstacle | |
McCowan | Microphone arrays: A tutorial | |
CN109490850B (zh) | 主瓣干扰下宽带阵列自适应波束形成方法 | |
CN112180329B (zh) | 一种基于阵元随机均匀分布球阵反卷积波束形成的汽车噪声源声成像方法 | |
Wang et al. | Combining superdirective beamforming and frequency-domain blind source separation for highly reverberant signals | |
CN104166120B (zh) | 一种声矢量圆阵稳健宽带mvdr方位估计方法 | |
CN115276746B (zh) | 基于交替方向惩罚的频率一致宽带波束形成方法及系统 | |
CN115038012B (zh) | 基于admm的麦克风阵列鲁棒频率不变波束形成方法 | |
CN107102298B (zh) | 基于迭代互耦校正的雷达协方差矩阵重构波束形成方法 | |
CN111693971B (zh) | 一种用于弱目标检测的宽波束干扰抑制方法 | |
Xu et al. | A modified differential beamforming and its application for DOA estimation of low frequency underwater signal | |
CN113889136A (zh) | 一种基于麦克风阵列的拾音方法、拾音装置及存储介质 | |
Yang et al. | On the design of flexible Kronecker product beamformers with linear microphone arrays | |
CN106842135B (zh) | 基于干扰加噪声协方差矩阵重构的自适应波束形成方法 | |
CN108880586A (zh) | 一种宽带弱信号增强方法与装置 | |
CN110457796B (zh) | 无规平面差分传声器阵列阵元位置的优化方法 | |
CN109541526B (zh) | 一种利用矩阵变换的圆环阵方位估计方法 | |
CN113673419B (zh) | 适用于强干扰环境的波束域快速稀疏贝叶斯方位估计方法 | |
CN111551892B (zh) | 一种稳健自适应波束形成方法及装置 | |
Frank et al. | Constant-beamwidth linearly constrained minimum variance beamformer | |
CN109669172B (zh) | 基于主瓣内强干扰抑制的弱目标方位估计方法 | |
CN116840785A (zh) | 一种变电站低频声源定位系统及方法 | |
CN115276745B (zh) | 一种基于交替方向惩罚的鲁棒波束形成方法及系统 | |
CN113673158B (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 |