CN113532617A - 一种长时波束相位统计特性的线谱检测方法 - Google Patents

一种长时波束相位统计特性的线谱检测方法 Download PDF

Info

Publication number
CN113532617A
CN113532617A CN202110786961.XA CN202110786961A CN113532617A CN 113532617 A CN113532617 A CN 113532617A CN 202110786961 A CN202110786961 A CN 202110786961A CN 113532617 A CN113532617 A CN 113532617A
Authority
CN
China
Prior art keywords
sequence
theta
spectrum
matrix
phase
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
CN202110786961.XA
Other languages
English (en)
Other versions
CN113532617B (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.)
National University of Defense Technology
Original Assignee
National University of Defense 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 National University of Defense Technology filed Critical National University of Defense Technology
Priority to CN202110786961.XA priority Critical patent/CN113532617B/zh
Publication of CN113532617A publication Critical patent/CN113532617A/zh
Application granted granted Critical
Publication of CN113532617B publication Critical patent/CN113532617B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G01MEASURING; TESTING
    • G01HMEASUREMENT OF MECHANICAL VIBRATIONS OR ULTRASONIC, SONIC OR INFRASONIC WAVES
    • G01H3/00Measuring characteristics of vibrations by using a detector in a fluid
    • G01H3/04Frequency
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01HMEASUREMENT OF MECHANICAL VIBRATIONS OR ULTRASONIC, SONIC OR INFRASONIC WAVES
    • G01H3/00Measuring characteristics of vibrations by using a detector in a fluid
    • G01H3/10Amplitude; Power
    • 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

Landscapes

  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Measurement Of Velocity Or Position Using Acoustic Or Ultrasonic Waves (AREA)

Abstract

本发明属于水声工程、海洋工程和声呐技术等领域,具体涉及一种长时波束相位统计特性的线谱检测方法,先通过水平阵列的波束形成获得水平各个方位各个频点输出的复波束序列,在对准目标方位的输出序列信噪比得到了空间增强;进而以第一快拍的复波束序列为参考,计算其它时刻的复波束序列与参考复波束序列的互谱,获得互谱序列。本发明的特点在于首先对水平阵列进行波束形成,在增强线谱信号的信噪比的同时一定程度上抑制了其它方向上的强干扰;同时,基于长时互谱序列的相位统计特性来区分线谱与各类宽带噪声,其检测门限易于确定且相对稳定。

Description

一种长时波束相位统计特性的线谱检测方法
技术领域
本发明属于水声工程、海洋工程和声呐技术等领域,具体涉及一种长时波束相位统计特性的线谱检测方法。
背景技术
舰船辐射噪声通常由连续谱和线谱叠加而成,其中线谱噪声由辅机机械的往复运动引起,其频率低,可远距离传播,且谱级比连续谱高10-25dB,因此线谱是水下目标远程探测和识别的重要特征。线谱检测主要采用功率谱估计、现代谱估计等方法,二者都以离散傅里叶变换(Discrete Fourier Transform,DFT)为基础采用假设检验类线谱检测方法,依据信号的统计特性构建检测器。低频线谱分析历程图(Low Frequency Analysis AndRecording,LOFAR)利用线谱信号具有时间连续性的特性,通过对一段接收数据连续进行短时傅里叶变换(Short-Time Fourier Transform,STFT)计算得到。LOFAR图中线谱检测的主要原理是利用人眼的视觉累积效应,对一定时间内连续的线谱信号进行频率的自动跟踪。
为了提升线谱检测能力,美国海军实验室的Wagstaff等人利用线谱信号的相位特征具有长时稳定性,而噪声信号的相位长时间观测具有随机性,发明了Wagstaff’sintegration silence processor(WISPR)及一系列派生处理技术,可对相位缓慢变化的线谱信号采用相位差分的方法进行相位补偿从而实现线谱信号的长时间积分,以获得足够的处理增益(参见:Exploiting phase fluctuations to improve temporal coherence,IEEE J.Ocean.Eng.,2004,29(2):498–510.以及The AWSUMfilter:A 20-dB gainfluctuation-based processor,IEEE J.Ocean.Eng.,1997,22(1):110–118)。解恺等人结合空域阵列的输出增益,通过对波束域输出信号的起伏相位差分对齐,有效提高目标线谱分量的检测信噪比增益的分辨能力和波达方位的估计能力(参见:基于目标辐射噪声的信号起伏检测算法研究,《电子与信息学报》,2013,35(4):844-851)。
但是现有检测方法存在如下问题:一方面,海洋环境背景噪声起伏或者是其它舰船干扰将使得线谱的信噪比发生起伏,当线谱的某些部分因信噪比起伏,导致信噪比过低而淹没在干扰之中时,该部分线谱将无法进行检测。虽然使用阵列通过波束形成技术可以抑制水平面上来自其他方向的干扰信号,在一定程度上提升了期望目标方向上波束域信号的信噪比,但是强干扰存在时仍然可能使期望目标线谱淹没在噪声中;另一方面,由于海洋波导的干涉效应,目标运动时接收器接收到的线谱功率也会发生起伏,起伏范围能达到十几甚至二十分贝,这导致线谱在能量高时能够被检测到,在能量低时无法被检测到。考虑到噪声信号相位随机抖动,长时间观测相位方差大,而线谱信号相位相对稳定,长时间观测相位方差非常小,即使存在由于目标或者接收器相对运动带来的多普勒频移,其线谱信号相位方差在一定的时间内仍然能保持一个较小值。
发明内容
本发明所要解决的技术问题是,针对现有线谱检测抗干扰能力不足且检测门限难以确定的缺点,提供了一种长时波束相位统计特性的线谱检测方法。该方法将目标方向上波束域输出频域信号的长时间相位二次统计特性(即相位方差)作为检测器的检验量,在获得了阵列增益提高输出信噪比的同时,能够根据线谱和噪声的相位统计特性来区分线谱与各类宽带噪声(包括泄露到线谱方向上的宽带干扰能量),具有检测门限易于确定且抗干扰能力强的优点。
本发明所采用的技术方案是:
一种长时波束相位统计特性的线谱检测方法,利用一条由M个水听器阵元构成的水听器阵列将接收到的声信号进行频域波束形成,并求解目标波束方位各个频率点上长时间的互谱相位方差作为门限值,提升线谱检测的性能。具体分为以下步骤:
S1,对水听器阵列M个阵元接收到的声信号矩阵沿着每一列分窗口进行傅里叶变换,获得频谱矩阵P(fi);具体步骤如下:
S1.1将水听器阵列M个阵元接收到的声信号表示为一个NT×M维的矩阵p:
Figure BDA0003159398090000021
其中矩阵p的每一行p1(t) p2(t) … pM(t)表示水听器阵列的M个阵元在t时刻接收到的声信号,t=1,2,…,NT;设时间采样率为fs,总观测时间为T,则每个阵元总采样点数为NT=floor[Tfs],这里floor表示向下取整操作。
S1.2任取一个不大于NT的正整数Nwin作为窗口长度,并任取一个小于Nwin的正整数Noverlap作为窗口重叠长度,将矩阵p的每一列数据划分为K=floor[(NT-Nwin)/(Nwin-Noverlap)+1]个窗口;
S1.3选取矩阵p每一列的第k个窗口,即第(k-1)(Nwin-Noverlap)+1行到第(k-1)(Nwin-Noverlap)+Nwin行的数据进行Nwin点的离散傅里叶变换,获得频谱矩阵P(fi),k=1,2,...,K:
Figure BDA0003159398090000031
为傅里叶变换后的频率轴向量,阵列最低工作频率fmin和最高工作频率fmax在fFFT向量中对应的编号分别为imin和imax,imin=floor[fminNwin/fs]+1,imax=floor[fmaxNwin/fs]+1;记录第m列第k个窗口傅里叶变换后的第i个频点fi处的值为Pm(fi,k),其中fi=(i-1)fs/Nwin,m=1,2,...,M,i=imin,imin+1,...,imax;依次遍历第m列的K个窗口,然后遍历所有列,获得一个K×M维的频谱矩阵P(fi):
Figure BDA0003159398090000032
其中矩阵P(fi)第k行的数据P(fi,k)称为第k个快拍的数据:P(fi,k)=[P1(fi,k),...,PM(fi,k)],k称为快拍数。
S2,对频谱矩阵P(fi)的第k行数据P(fi,k)进行波束形成,获得波束形成后的复波束序列b(fi,θ,k),然后计算水听器阵列工作频段内宽带波束功率输出结果B(θ,k)及其极大值对应的可能目标方位角;具体步骤如下:
S2.1对P(fi,k)进行波束形成,获得复波束序列
b(fi,θ,k)=wH(fi,θ)PT(fi,k)
式中上标H表示矩阵的复共轭转置,上标T表示矩阵的转置;w(fi,θ)为远场平面波模型下水听器阵列的导向向量;θ为水平方位搜索角,角度搜索间隔为Δθ,其最小值和最大值分别为θmin和θmax
S2.2根据S2.1获得的复波束序列b(fi,θ,k),求解水听器阵列工作频段内宽带波束功率输出结果:
Figure BDA0003159398090000033
S2.3记录所有使B(θ,k)达到极大值时对应的方位角θsk,s=1,2,...,S,S表示B(θ,k)极大值的个数,即第k个快拍下可能的目标个数。
S3,计算在第s个极大值方位角θsk处、频点为fi时,复波束序列b(fi,θ,k)在第k个快拍的序列值与第1个快拍的序列值的复互谱值C(fi,s,k)及其相位值φ(fi,s,k):
C(fi,s,k)=b(fisk,k)b*(fis1,1)
φ(fi,s,k)=arg{C(fi,s,k)}
其中上标*表示向量的复共轭,arg{·}表示对复数取相位的操作。
S4,重复S3直至遍历所有K个快拍数,获得长度为K的长时统计复互谱相位序列φ(fi,s),然后计算该序列的方差δφ(fi,s)。
长时统计复互谱相位序列可以表示为
φ(fi,s)=arg{[C(fi,s,1),C(fi,s,2),...,C(fi,s,K)]T}
其方差为
Figure BDA0003159398090000041
其中
Figure BDA0003159398090000042
为第i个频点的平均相位。
S5,将第imin至imax共计imax-imin+1个频点处所有相位方差δφ(fi,s)的倒数的平均值
Figure BDA0003159398090000043
作为检测门限值,然后将不同频点处δφ(fi,s)的倒数逐个与
Figure BDA0003159398090000044
比较:超出门限值的视为线谱,低于门限的值视为随机噪声。其中门限值
Figure BDA0003159398090000045
可表示为
Figure BDA0003159398090000046
S6,在第s个极大值方位角θsk处,重复S3至S5直至遍历所有imax-imin+1个频点,得到第s个可能目标的线谱检测结果。
S7,重复S6直至遍历所有s个可能目标,即可完成线谱检测。
本发明基于以下原理:线谱信号的相位具有长时稳定性,而噪声信号相位的长时间序列是随机的,本发明先通过水平阵列的波束形成获得水平各个方位各个频点输出的复波束序列,在对准目标方位的输出序列信噪比得到了空间增强。进而以第一快拍的复波束序列为参考,计算其它时刻的复波束序列与参考复波束序列的互谱,获得互谱序列。对于线谱信号,该互谱序列的相位反映了目标波束方向上线谱信号的相位长时稳定性。对于噪声或者其它宽带信号,该互谱序列的相位长时间统计结果是随机的。本发明的特点在于首先对水平阵列进行波束形成,在增强线谱信号的信噪比的同时一定程度上抑制了其它方向上的强干扰;同时,基于长时互谱序列的相位统计特性来区分线谱与各类宽带噪声(包括泄露到线谱方向上的宽带干扰能量),其检测门限易于确定且相对稳定。
与现有技术相比,本发明所具有的有益效果为:
本发明通过对水平阵列进行波束形成并计算互谱相位方差获得了线谱检测的门限值。与传统的方法相比,本发明的门限获取方法具有自适应性,且计算量小,稳健性好,便于工程实现;相比传统的功率谱法能够获得更高的线谱信噪比。
附图说明
图1是本发明实施例的总体流程图;
图2是阵列波束形成方位历程结果图;
图3是对准期望目标方向波束输出的LOFAR结果图;
图4是使用目标波束方向上的LOFAR输出检测结果(a)和本发明实施例的波束互谱输出的相位统计特性检测结果(b)比对图,其中五角星表示声源发射的线谱信号。
具体实施方式
下面结合附图和某次海试数据处理的具体实施方式,对本发明进行详细阐述。本发明提出的一种长时波束相位统计特性的线谱检测方法,通过对波束输出互谱相位进行统计分析,当检测到长时相位异常判别为线谱。附图1给出本发明实施例的总体流程图,具体将结合某次海上试验得到的数据进行详细说明。
在此次海上试验中,采用水听器个数M=28的不等间距阵列,其阵列形状接近线列阵。以阵列第1个阵元为原点,建立xoy直角坐标系,其第m个阵元的坐标为xm=[xm,ym]T。水平方位搜索角θ定义为与x正轴的夹角,角度搜索范围从0度到359度,角度搜索间隔Δθ为1度。
S1,对水听器阵列M个阵元接收到的声信号矩阵沿着每一列分窗口进行傅里叶变换,获得频谱矩阵P(fi);具体步骤如下:
S1.1将水听器阵列M个阵元接收到的声信号表示为一个NT×M维的矩阵p:
Figure BDA0003159398090000051
其中矩阵p的每一行p1(t) p2(t) … pM(t)表示水听器阵列的M个阵元在t时刻接收到的声信号,t=1,2,…,NT;设时间采样率为fs,总观测时间为T,则每个阵元总采样点数为NT=floor[Tfs],这里floor表示向下取整操作。
S1.2任取一个不大于NT的正整数Nwin作为窗口长度,并任取一个小于Nwin的正整数Noverlap作为窗口重叠长度,将矩阵p的每一列数据划分为K=floor[(NT-Nwin)/(Nwin-Noverlap)+1]个窗口;
S1.3选取矩阵p每一列的第k(k=1,2,...,K)个窗口,即第(k-1)(Nwin-Noverlap)+1行到
第(k-1)(Nwin-Noverlap)+Nwin行的数据进行Nwin点的离散傅里叶变换,获得频谱矩阵P(fi):
Figure BDA0003159398090000061
为傅里叶变换后的频率轴向量,阵列最低工作频率fmin和最高工作频率fmax在fFFT向量中对应的编号分别为imin和imax,那么imin=floor[fminNwin/fs]+1,imax=floor[fmaxNwin/fs]+1;记录第m列(m=1,2,...,M)第k个窗口傅里叶变换后的第i个频点fi处的值为Pm(fi,k),其中fi=(i-1)fs/Nwin(i=imin,imin+1,...,imax);依次遍历第m列的K个窗口,然后遍历所有列,可获得一个K×M维的频谱矩阵P(fi),即
Figure BDA0003159398090000062
其中矩阵P(fi)第k行的数据称为第k个快拍的数据,k称为快拍数。
本实施例中,时间采样率fs=3276.8Hz,选取总观测时间为T=5min,每个水听器总采样点数为NT=983040。进行傅里叶变换时采用的时间窗口点数为Nwin=65536个点,即每个时间窗口的观测时长为20s。每次选取数据时与上一次数据的重叠点数Noverlap=45875,因此总采样时间时间T内,一共可以划分为K=47个窗口。因此,矩阵P(fi)的维数为47×28。由于时间窗口点数为65536个点,故傅里叶变换后频点数为65536个。阵列工作频段为40Hz到410Hz,在频率轴向量fFFT上对应的频点编号分别为imin=802和imax=8202,因此总的频率点数为imax-imin+1=7401个。
S2,对频谱矩阵P(fi)的第k行数据进行波束形成,获得波束形成后的复波束序列b(fi,θ,k)和其工作频段内的宽带波束功率输出结果B(θ,k)及其极大值对应的可能目标方位角;具体步骤如下:
S2.1将P(fi,k)=[P1(fi,k),...,PM(fi,k)]记为频谱矩阵P(fi)的第k行数据,即水听器阵列M个阵元第k个快拍的声信号经过傅里叶变换后得到第i个频点的频谱信号。对P(fi,k)进行波束形成,获得复波束序列
b(fi,θ,k)=wH(fi,θ)PT(fi,k)
式中上标H表示矩阵的复共轭转置,T表示矩阵的转置;θ为水平方位搜索角,角度搜索间隔为Δθ,其最小值和最大值分别为θmin和θmax,因此θ=[θminmin+Δθ,...,θmax]T;不失一般性地,设θmin=-180°,θmax=180°,在实际操作中,可根据阵列对目标的有效探测方位角范围这一先验信息缩小水平方位搜索角的范围,例如阵列由于存在左右模糊只能进行[-90°,90°]入射角范围的波束形成,可以将θmin和θmax分别设置为-90°和90°;w(fi,θ)为远场模型下平面波对应的水听器阵列导向向量。通常,在远场模型下平面波波束形成对应的水听器阵列的导向向量w(fi,θ)可表示为
Figure BDA0003159398090000071
为方便分析,后续内容将以上述导向向量为例。类似地,w(fi,θ)也可以为非等间距阵列或者其它几何构型的阵列在远场模型下平面波波束形成对应的水听器阵列的导向向量,不同类型的水听器阵列导向向量的表达式参见文献:R.O.Nielsen.Sonar SignalProcessing.Artech House,London,1991。
S2.2根据S2.1获得的复波束序列,进一步求解其阵列工作频段内宽带功率输出结果
Figure BDA0003159398090000072
S2.3记录所有使B(θ,k)达到极大值对应的方位角θsk(s=1,2,...,S),这里S表示B(θ,k)极大值的个数,即第k个快拍下可能的目标个数。优选地,本领域技术人员可根据期望目标方位的先验信息(例如通过船舶自动识别系统(AIS)获得期望目标方位的先验信息)减少可能的目标个数(如假设该步骤中记录到的极大值有3个,分别对应方位角θ1k、θ2k、θ3k分别为79度、80度和100度,而根据实验过程的先验信息表明可能的目标只会存在于80度附近,当方位角排除门限设为5度时,可排除100度方位角对应的极大值,从而将可能的目标个数减少为2个)。
本实施例中,水平方位角空间扫描角度为1到180度,角度间隔为1度。因此,获得了不同频率点fi和不同快拍k情况下,180个角度方向上的复波束序列输出b(fi,θ,k)。按照波束序列的功率输出B(θ,k),得到了附图2所示的目标方位历程图。附图2给出了5分钟内B(θ,k)的输出结果。根据试验条件中的先验信息,期望目标位于80度附近。从附图2可以看到,在5分钟内,强干扰方位在8度左右,期望目标角度位于78度基本不变,因此θ11=...=θ1K=78度。根据附图1,干扰能量高于目标。如附图3所示,期望目标方向上输出波束信号特征,即LOFAR分析结果|b(fi1k,k)|2。比对附图2和附图3可知,虽然干扰方向和目标方向差异很大,但是强干扰的能量仍然通过旁瓣泄露到了目标方向,导致此目标方位上LOFAR图的低频部分存在很强的舰船干涉能量,而且很多线谱能量不稳定,无法从LOFAR图中观察到。
S3,计算在第s个极大值方位角θsk处、频点为fi时,复波束序列在第k个快拍的序列值与第1个快拍的序列值的复互谱值C(fi,s,k)及其相位值φ(fi,s,k):
C(fi,s,k)=b(fisk,k)b*(fis1,1)
φ(fi,s,k)=arg{C(fi,s,k)}
其中上标*表示向量的复共轭,arg{}表示对复数取相位的操作。
本实施例中,由于数据点数被划分为K=47个窗口,总的频率点数为imax-imin+1=7401个,因此共获得47个复互谱相位值φ(fi,k),即每47个快拍进行后续的线谱检测。因此对阵列工作频带范围内的7401个频点,获取的互谱相位φ(fi,k)均为47个。
S4,重复S3直至遍历所有K个快拍数,获得长度为K的长时统计复互谱相位序列,然后计算该序列的方差δφ(fi,s,k)。
长时统计复互谱相位序列可以表示为
φ(fi,s,k)=arg{[C(fi,s,1),C(fi,s,2),...,C(fi,s,K)]T}
其方差为
Figure BDA0003159398090000091
其中
Figure BDA0003159398090000092
为第i个频点的平均相位。
S5,将第imin至imax共计imax-imin+1个频点处所有相位方差δφ(fi,s)的倒数的平均值
Figure BDA0003159398090000093
作为检测门限值,然后将不同频点处δφ(fi,s)的倒数逐个与
Figure BDA0003159398090000094
比较:超出门限值的视为线谱,低于门限的值视为随机噪声。其中门限值
Figure BDA0003159398090000095
可表示为
Figure BDA0003159398090000096
S6,在第s个极大值方位角θsk处,重复S3至S5直至遍历所有imax-imin+1个频点,得到第s个可能目标的线谱检测结果。
S7,重复S6直至遍历所有s个可能目标,即可完成线谱检测。
本实施例中,对阵列工作频带范围内的每个频点,获得期望目标的互谱相位φ(fi,1,k)的均值
Figure BDA0003159398090000097
和方差δφ(fi,1)。附图4(b)给出了所有频点的相位方差δφ(fi,1)。为了验证本发明的有效性,图4给出了传统的功率谱法(附图4(a))和本发明的结果(附图4(b))比对。传统的功率谱法通过将附图3中所示的期望目标方向上输出波束信号特征进行多个快拍的平均,即
Figure BDA0003159398090000098
在本实例中存在宽带强干扰目标且信干比低,此时传统的功率谱法难以通过信号能量检测到线谱。即使对附图4(a)通过去趋势平滑功率谱的趋势项,也无法提高线谱的信噪比,依然难以检测到信号。本发明利用线谱信号的相位稳定性的统计规律与宽带噪声随机相位的统计规律,将线谱信号与噪声分离,而不是采用信号能量特征,具有更好的检测效果。同时,不同噪声频点的相位方差非常接近,即噪声相位方差谱非常平坦,这为门限确定提供了便利,便于系统实现。
最后应说明的是:以上所述仅为本发明的较佳实施例而已,并非用于用于限定本发明的保护范围,尽管参照前述实施例对本发明进行了详细的说明,对于本领域的技术人员来说,其依然可以对前述各实施例所记载的技术方案进行修改,或者对其中部分技术特征进行等同替换。凡在本发明的精神和原则之内,所作的任何修改、等同替换、改进等,均应包含在本发明的保护范围之内。

Claims (3)

1.一种长时波束相位统计特性的线谱检测方法,其特征在于:该方法利用一条由M个水听器阵元构成的水听器阵列将接收到的声信号进行频域波束形成,并求解目标波束方位各个频率点上长时间的互谱相位方差作为门限值,提升线谱检测的性能,具体分为以下步骤:
S1,对水听器阵列M个阵元接收到的声信号矩阵沿着每一列分窗口进行傅里叶变换,获得频谱矩阵P(fi);具体步骤如下:
S1.1将水听器阵列M个阵元接收到的声信号表示为一个NT×M维的矩阵p:
Figure FDA0003159398080000011
其中矩阵p的每一行p1(t) p2(t) … pM(t)表示水听器阵列的M个阵元在t时刻接收到的声信号,t=1,2,…,NT;设时间采样率为fs,总观测时间为T,则每个阵元总采样点数为NT=floor[Tfs],这里floor表示向下取整操作;
S1.2任取一个不大于NT的正整数Nwin作为窗口长度,并任取一个小于Nwin的正整数Noverlap作为窗口重叠长度,将矩阵p的每一列数据划分为K=floor[(NT-Nwin)/(Nwin-Noverlap)+1]个窗口;
S1.3选取矩阵p每一列的第k个窗口,即第(k-1)(Nwin-Noverlap)+1行到第(k-1)(Nwin-Noverlap)+Nwin行的数据进行Nwin点的离散傅里叶变换,获得频谱矩阵P(fi),k=1,2,...,K:
Figure FDA0003159398080000012
为傅里叶变换后的频率轴向量,阵列最低工作频率fmin和最高工作频率fmax在fFFT向量中对应的编号分别为imin和imax,imin=floor[fminNwin/fs]+1,imax=floor[fmaxNwin/fs]+1;记录第m列第k个窗口傅里叶变换后的第i个频点fi处的值为Pm(fi,k),其中fi=(i-1)fs/Nwin,m=1,2,...,M,i=imin,imin+1,...,imax;依次遍历第m列的K个窗口,然后遍历所有列,获得一个K×M维的频谱矩阵P(fi):
Figure FDA0003159398080000021
其中矩阵P(fi)第k行的数据P(fi,k)称为第k个快拍的数据:P(fi,k)=[P1(fi,k),...,PM(fi,k)],k称为快拍数;
S2,对频谱矩阵P(fi)的第k行数据P(fi,k)进行波束形成,获得波束形成后的复波束序列b(fi,θ,k),然后计算水听器阵列工作频段内宽带波束功率输出结果B(θ,k)及其极大值对应的可能目标方位角;具体步骤如下:
S2.1对P(fi,k)进行波束形成,获得复波束序列
b(fi,θ,k)=wH(fi,θ)PT(fi,k)
式中上标H表示矩阵的复共轭转置,上标T表示矩阵的转置;w(fi,θ)为远场平面波模型下水听器阵列的导向向量;θ为水平方位搜索角,角度搜索间隔为Δθ,其最小值和最大值分别为θmin和θmax
S2.2根据S2.1获得的复波束序列b(fi,θ,k),求解水听器阵列工作频段内宽带波束功率输出结果:
Figure FDA0003159398080000022
S2.3记录所有使B(θ,k)达到极大值时对应的方位角θsk,s=1,2,...,S,S表示B(θ,k)极大值的个数,即第k个快拍下可能的目标个数;
S3,计算在第s个极大值方位角θsk处、频点为fi时,复波束序列b(fi,θ,k)在第k个快拍的序列值与第1个快拍的序列值的复互谱值C(fi,s,k)及其相位值φ(fi,s,k):
C(fi,s,k)=b(fisk,k)b*(fis1,1)
φ(fi,s,k)=arg{C(fi,s,k)}
其中上标*表示向量的复共轭,arg{·}表示对复数取相位的操作;
S4,重复S3直至遍历所有K个快拍数,获得长度为K的长时统计复互谱相位序列φ(fi,s),然后计算该序列的方差δφ(fi,s);
长时统计复互谱相位序列可以表示为
φ(fi,s)=arg{[C(fi,s,1),C(fi,s,2),...,C(fi,s,K)]T}
其方差为
Figure FDA0003159398080000031
其中
Figure FDA0003159398080000032
为第i个频点的平均相位;
S5,将第imin至imax共计imax-imin+1个频点处所有相位方差δφ(fi,s)的倒数的平均值
Figure FDA0003159398080000033
作为检测门限值,然后将不同频点处δφ(fi,s)的倒数逐个与
Figure FDA0003159398080000034
比较:超出门限值的视为线谱,低于门限的值视为随机噪声;其中门限值
Figure FDA0003159398080000035
可表示为
Figure FDA0003159398080000036
S6,在第s个极大值方位角θsk处,重复S3至S5直至遍历所有imax-imin+1个频点,得到第s个可能目标的线谱检测结果;
S7,重复S6直至遍历所有s个可能目标,即可完成线谱检测。
2.一种根据权利要求1所述长时波束相位统计特性的线谱检测方法,其特征在于:S2.1中,远场模型下平面波波束形成对应的水听器阵列的导向向量w(fi,θ)可表示为
Figure FDA0003159398080000037
3.一种根据权利要求1所述长时波束相位统计特性的线谱检测方法,其特征在于:w(fi,θ)也可以为非等间距阵列或者其它几何构型的阵列在远场模型下平面波波束形成对应的水听器阵列的导向向量。
CN202110786961.XA 2021-07-13 2021-07-13 一种长时波束相位统计特性的线谱检测方法 Active CN113532617B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202110786961.XA CN113532617B (zh) 2021-07-13 2021-07-13 一种长时波束相位统计特性的线谱检测方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202110786961.XA CN113532617B (zh) 2021-07-13 2021-07-13 一种长时波束相位统计特性的线谱检测方法

Publications (2)

Publication Number Publication Date
CN113532617A true CN113532617A (zh) 2021-10-22
CN113532617B CN113532617B (zh) 2023-11-03

Family

ID=78127579

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202110786961.XA Active CN113532617B (zh) 2021-07-13 2021-07-13 一种长时波束相位统计特性的线谱检测方法

Country Status (1)

Country Link
CN (1) CN113532617B (zh)

Citations (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JPH1123700A (ja) * 1997-07-01 1999-01-29 Tech Res & Dev Inst Of Japan Def Agency ノーマルモード固有値推定方法
US5914912A (en) * 1997-11-28 1999-06-22 United States Of America Sonar array post processor
US20060133211A1 (en) * 2004-12-17 2006-06-22 Yang Tsih C Method and apparatus for acoustic source tracking using a horizontal line array
CN104052535A (zh) * 2014-06-23 2014-09-17 东南大学 基于空分多址与干扰抑制的毫米波大规模mimo系统多用户传输方法
CN105301580A (zh) * 2015-10-30 2016-02-03 哈尔滨工程大学 一种基于分裂阵互谱相位差方差加权的被动探测方法
CN107179535A (zh) * 2017-06-01 2017-09-19 东南大学 一种基于畸变拖曳阵的保真增强波束形成的方法
CN110138422A (zh) * 2019-04-16 2019-08-16 电子科技大学 基于稀疏编码和无相位解码的毫米波通信快速波束对准方法
CN111665489A (zh) * 2019-03-08 2020-09-15 中国科学院声学研究所 一种基于目标特性的线谱提取方法

Patent Citations (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JPH1123700A (ja) * 1997-07-01 1999-01-29 Tech Res & Dev Inst Of Japan Def Agency ノーマルモード固有値推定方法
US5914912A (en) * 1997-11-28 1999-06-22 United States Of America Sonar array post processor
US20060133211A1 (en) * 2004-12-17 2006-06-22 Yang Tsih C Method and apparatus for acoustic source tracking using a horizontal line array
CN104052535A (zh) * 2014-06-23 2014-09-17 东南大学 基于空分多址与干扰抑制的毫米波大规模mimo系统多用户传输方法
CN105301580A (zh) * 2015-10-30 2016-02-03 哈尔滨工程大学 一种基于分裂阵互谱相位差方差加权的被动探测方法
CN107179535A (zh) * 2017-06-01 2017-09-19 东南大学 一种基于畸变拖曳阵的保真增强波束形成的方法
CN111665489A (zh) * 2019-03-08 2020-09-15 中国科学院声学研究所 一种基于目标特性的线谱提取方法
CN110138422A (zh) * 2019-04-16 2019-08-16 电子科技大学 基于稀疏编码和无相位解码的毫米波通信快速波束对准方法

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
朱秀珍;乔长城;杨清海;: "基于矢量水听器瞬时相位差方差加权的低频线谱检测算法", 电声技术, no. 01 *

Also Published As

Publication number Publication date
CN113532617B (zh) 2023-11-03

Similar Documents

Publication Publication Date Title
CN112269164B (zh) 深海可靠声路径下基于干涉结构匹配处理弱目标定位方法
CN110196414B (zh) 一种基于补偿天线方向图误差的天线波束指向方法
CN109459744B (zh) 一种实现多干扰抑制的稳健自适应波束形成方法
CN108872970B (zh) 适用于一般等间距稀疏阵列单频信号波束形成的栅瓣判别方法
CN109541548B (zh) 一种基于匹配场的空气声呐定位方法
JP2011158471A (ja) 時空間適応処理システムにおいてターゲットを検出するための方法
CN111693971B (zh) 一种用于弱目标检测的宽波束干扰抑制方法
CN110687528B (zh) 自适应波束形成器生成方法及系统
CN110161489A (zh) 一种基于伪框架的强弱信号测向方法
CN112881975B (zh) 基于子阵特征矩阵联合对角化的单脉冲和差波束测角方法
US7319640B1 (en) Noise suppression system
CN110261814B (zh) 基于空间谱重构和导向矢量直接估计的波束形成方法
CN113532617A (zh) 一种长时波束相位统计特性的线谱检测方法
Wang et al. Snapshot-deficient active target localization in beam-time domain using multi-frequency expectation-maximization algorithm
CN109814065B (zh) 基于相位因子加权的波束形成方法
CN111665469A (zh) 一种基于空间时频分布的水下多径信号参数估计方法
CN110850421A (zh) 基于混响对称谱的空时自适应处理的水下目标检测方法
CN115561764A (zh) 一种基于单矢量水听器的运动目标深度估计方法
CN111414580B (zh) 一种低信混比条件下的混响抑制方法
CN110632579B (zh) 一种利用子阵波束域特征的迭代波束形成方法
CN114647931A (zh) 一种基于期望信号消除和空间谱估计的稳健波束形成方法
Ma et al. Passive direction-of-arrival estimation under high noise and strong interference condition for volumetric array
CN110007296A (zh) 一种基于引导信号修正的时域干扰抵消方法
CN112114287B (zh) 一种方位观测数据的野值实时剔除方法
Zhang et al. Fast Estimation of Direction of Arrival for Towed Array Based on Sparse Bayesian Learning

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