专利名称：一种快速宽带频域扩展拖曳阵波束形成方法
技术领域：
本发明属于声纳数字信号处理领域，特别涉及一种快速宽带频域扩展拖曳阵波束形成方法。

背景技术：
被动合成孔径(PASA)是一种被动接收目标噪声、通过阵列的运动和信号处理方法人工增加小孔径阵列长度的技术，可以使人们依靠短阵的机动获得长阵的各种好处。如今，被动声纳已经广泛使用在众多海军系统中，如拖曳阵、舷侧阵、鱼类声纳以及海底线列阵等，主要用来定位、检测、分类，或者对目标进行跟踪。近几年来，随着人们对被动合成孔径的关注程度逐渐提高，其性能被不断加强。
现有的定位主要通过波束形成来实现，它利用目标信号到达不同阵元传感器的时间差，通过相位补偿融合了每个传感器的接收信号，以提高输出信号的信噪比来估计平面波信号的波达方向。实际上，目标辐射信号一般都是宽带信号，常规做法首先将阵列信号进行时域FFT变换，然后针对每个频点的响应做单频点波束形成，最后将所有频点的波束结果融合在一起，得到最终波束输出。
对于常规波束形成，输出信噪比、方位分辨率以及声纳作用距离随着阵长度的增加而增大。但实际上，受平台机动以及各种因素限制，线列阵不能设计的太长，从而制约了声纳的方位分辨率和作用距离。
随着计算性能的增加，通过处理阵列在空间和时间上接收的信号，根据已有小孔径阵列虚拟出大孔径阵列已经成为可能，这就是现在被广泛研究的被动合成孔径技术。阵列做匀速直线运动时，可以利用阵列在不同时间间隔上接收到信号存在的相位差，进行信号相关融合，从而生成一个大的虚拟孔径。实际上，由于积分时间内存在介质扰动，而且阵列在运动路径内也会出现不规则扰动，导致一定时间间隔内阵元接收信号的相位差随之出现不确定的扰动，因此每进行一次被动孔径合成都需要对相位变化进行补偿。
常用的被动合成孔径技术有两种扩展拖曳阵技术(ETAM)，如文献1“S.Stergiopoulos，E.H.Sullivan，Extended towed array processing by an overlap correlator，J.Acoust.Soc.Am.，vol.86，pp.158-171，1989”；快速傅里叶变换被动合成孔径技术(FFTSA)，如文献2“Stergiopoulos and Urban，A new passive synthetic aperturetechnique for towed arrays，IEEE Journal of Oceanic Engineering，Vol 17，No.1，January1992”。
ETAM使用重叠相关方法来估计快拍之间由介质和路径扰动引起的相位扰动。文献3“S.Stergiopoulos，Optimum bearing resolution for a moving towed array andextension of its physical aperture，J.Acoust.Soc.Am.，vol.87，pp.2128-2140，1990”评估了ETAM算法的性能，它和最大似然估计器类似性能相当，当水听器接收信号的信噪比大于0dB时，ETAM算法能够提供和实际阵列相比拟的波束估计。
FFTSA算法直接将一定时间间隔内两个快拍信号的相位差作为常数，而没有考虑由于介质或运动偏差导致的相位扰动。其好处是简化了计算，使得孔径合成过程可以直接用FFT算法代替进行，大大提高了计算速度；但是，速度的提高却是以降低性能为代价换取的。
文献4“R.Rajagopal，P.Ramakrishna Rao，Performance comparison of PASAbeamforming algorithms，International Symposium on Signal Processing and itsApplications，ISSPA，Gold Coast，Australia，1996”表明，比较ETAM、FFTSA和最大似然估计三种算法，无论在提高阵增益、检测概率(对于给定的虚警率)方面，还是针对不同输入信噪比时空时相关的稳定性，ETAM算法都表现出更优异的性能。另外，对于一个具有稳定动态声源谱的移动目标，仿真数据表明ETAM算法比FFTSA具有更好的稳定性。
总体而言，ETAM算法比FFTASA算法具有无可比拟的性能优势，但计算速度较慢，可行性不高。实际应用中需要一种能够针对宽带目标信号快速合成孔径、并进行波束形成以获得目标方位的实时处理方法。


发明内容
本发明目的在于，为克服现有的ETAM算法计算速度较慢，可行性不高且实际应用中需要一种能够针对宽带目标信号快速合成孔径、并进行波束形成以获得更精确目标方位的实时处理方法，从而提出一种快速宽带频域扩展拖曳阵波束形成方法。
本发明所述的一种快速宽带频域扩展拖曳阵波束形成方法，将ETAM被动合成孔径技术和宽带频域波束形成方法结合在一起，形成了一种新的基于ETAM被动合成孔径技术的宽带频域快速波束形成方法。该方法比常规波束形成方法具有更高的方位分辨率，对微弱目标的检测能力也更优异，同时速度较快，便于实时处理。
一种线列阵声纳装置，可以是一条拖曳阵，也可以是舷侧阵，所述线列阵由多个水听器组成。这里设实际阵元数目N，阵元间距d，阵做匀速直线运动，速度v；目标入射方向θ，阵元接收信号表示为x(t)；快拍长度为L。
为了实现上述目的，本发明提出一种快速宽带频域扩展拖曳阵波束形成方法，该方法将ETAM被动合成孔径技术和宽带频域波束形成方法结合在一起，可用于对微弱目标的检测和定位，同时速度较快，便于实时处理；所述的方法包含如下步骤 1)用线阵接收空间信号，得到M个阵元的时域信号；取重叠阵元数目为N，其中N大于零且小于阵元数目，计算两次快拍数据的时间间隔阵τ，即阵运动N个阵元距离需要的时间，如图1所示； 2)对第k快拍和第k+1快拍数据在时间域上做快速傅里叶变换； 3)确定目标辐射信号的频带范围[fmin，fmax]，根据第k快拍和第k+1快拍的重叠阵元计算相移参数ψk+1(fi) k表示快拍序号，n表示阵元序号，fi表示频率，ρm(fi)表示第n阵元频率fi的补偿因子 4)对第k+1快拍、第M-N+1～M阵元信号的不同频率分量分别进行相位补偿，将其虚拟成第k快拍的第M+1～M+N阵元；依次类推，这N个阵元信号可以虚拟成为第1快拍的第M+N*(k-1)+1～M+N*K阵元； 5)重复以上步骤，进行K次的孔径合成，可以虚拟得到MM＝M+N*K个阵元接收的信号在[fmin，fmax]范围内的频域响应； 6)对MM个阵元信号在空间域上补零，再对其在空间域上做FFT变换，然后将获得数据的零频分量移至谱中心； 7)在处理带宽内的每一个频点上依据频率——波束网格对其进行校正，得到不同频点的波束输出； 8)合成不同频率的波束形成结果，得到MM个阵元在不同方位上的功率输出，从而实现目标的精确定位。
上述技术方案，步骤1)所述的重叠阵元占阵元总数的一半，能取得更好的方向分辨率和技术效果。
作为本发明的一个改进，所述的空间域FFT运算结果用频率——波束网格对其进行校正得到单频信号的波束输出； 作为本发明的又一改进，步骤6)所述的空间域采用补零法，所述的补零方法具体采用直接在原数据后直接补零，或在原数据中插值补零。
在上述技术方案的基础上，进一步地，所述步骤6中，为了提高运算速度，采用FFT代替来常规波束形成。如公式(7)所示，频率fi的波束形成就是对所有阵元信号的单频响应进行相位补偿后做叠加运算，得到不同方位的波束输出。此过程与FFT运算非常相似，公式(8)表示对阵元信号单频响应做FFT运算。
对比公式(7)、(8)可以看出，波束方位θs和FFT变换域k存在一定对应关系 即频率——波束网格，这里k是FFT的运算序号，fi表示频率，d是阵元间距，θs是扫描方位，取值范围0～180度，MM是阵元数目。也就是说波束输出L(fi，θs)和阵列信号的空间域FFT运算结果Y(fi，k)存在对应关系，只要按照频率——波束网格对阵列信号的空间域FFT运算结果进行校正就能得到单频信号的波束输出。
空间域上补零的目的是保证阵元数目远大于波束数目，如公式(9)所示，当阵元数目MM足够大时，θs能够平均分布在不同的k值上，这样根据空间域FFT运算得到的波束输出也就越精确。空间域补零可以有多种选择。既可以直接在原数据后直接补零，也可以在原数据中插值补零，只要补零后的数据长度满足要求即可。
在上述技术方案的基础上，进一步地，所述步骤7中的频率——波束网格校正是预先产生的，根据公式(9)，计算方法如下所示，floor表示取整运算。
采用FFT实现频域波束形成需要依据频率——波数网格进行校正，该问题在文献“Brian Maranda，Efficient digital beamforming in the frequency domain，1989，J.Acoustical Society of America”中有详细的阐述。对本领域技术人员来说，理解和实现该方法是可以胜任的。
本发明的优点在于 (1)将ETAM被动合成孔径技术和宽带频域波束形成方法结合在一起，形成了一种新的基于ETAM被动合成孔径技术的宽带频域快速波束形成方法，比常规波束形成方法具有更高的方位分辨率，对微弱目标的检测能力也更优异； (2)ETAM孔径合成和波束形成运算都在频域实现，采用快速傅立叶变换，将宽带信号分解为多个频点信号，只对信号频带范围内的每个频点进行处理，速度快，实时性好，算法稳健性高； (3)海试数据证明有效。



图1是现有的ETAM方法被动合成孔径原理示意图； 图1-1本发明实施例中的ETAM方法被动合成孔径原理示意图； 图2是现有技术中由舰船拖曳的常规线列阵示意图； 图3是本发明的算法流程图； 图4是采用本发明算法与未采用本发明算法处理得到目标波束输出对比； 图5采用本发明处理结果和对同样长度实孔径阵列的波束处理结果对比； 图6是对采用常规实孔径波束形成处理海试数据得到的目标方位历程； 图7是对采用本发明处理海试数据得到的目标方位历程。

具体实施例方式 下面结合附图对本发明进行进一步说明。
本发明的基本构思是常规阵列由于孔径太短，输出信噪比较小，导致方位分辨率不高，且作用距离较短。当阵列做匀速直线运动时，ETAM被动合成孔径方法能够根据相隔一定延时的两个快拍信号之间的重叠阵元，计算由于传播介质和运动误差引起的相位扰动，然后对阵元信号进行相位补偿从而虚拟出大孔径阵元，提高输出信噪比，获得更高的方位分辨率和更远的作用距离。但ETAM算法比较耗时，难以实时性处理。本发明直接在频域对宽带目标信号进行ETAM孔径合成，然后利用宽带频域快速波束形成算法对合成的阵列数据做波束形成，得到不同方位的功率输出，从而快速精确得到远距离目标的方位。
本发明提出的一种新的基于ETAM被动合成孔径技术的宽带频域快速波束形成方法，包括如下步骤 步骤1用线阵接收空间信号，得到M个阵元的时域信号；取重叠阵元数目为

计算两次快拍数据的时间间隔阵τ，即阵运动

个阵元距离需要的时间，如图1-1。
步骤2对第k、k+1快拍数据在时间域上做快速傅里叶变换，如第k快拍第m阵元的时域信号变换过程如下 步骤3确定目标辐射信号的频带范围[fmin，fmax]，根据第k快拍和第k+1快拍的重叠阵元计算相移参数ψk+1(fi) k表示快拍序号，m表示阵元序号，fi表示频率，ρm(fi)表示第m阵元频率fi的补偿因子 步骤4对第k+1快拍、第

阵元信号的不同频率分量分别进行相位补偿，将其虚拟成第k快拍的第

阵元。依次类推，这

个阵元信号可以虚拟成为第1快拍的第

阵元。也就是说每进行一次孔径合成，可以虚拟出

个阵元。
步骤5重复以上步骤，进行K次的孔径合成，可以虚拟得到

个阵元接收的信号在[fmin，fmax]范围内的频域响应。
到这一步为止，合成孔径已经结束，下面进行快速宽带频域波束形成。
步骤6对MM个阵元信号在空间域上补零，再对其在空间域上做快速傅立叶变换，然后将获得数据的零频分量移至谱中心。
步骤7在处理带宽内的每一个频点上依据频率——波束网格对其进行校正，得到不同频点的波束输出。
步骤8合成不同频率的波束形成结果，得到MM个阵元在不同方位上的功率输出。
在上述技术方案的基础上，进一步地，所述步骤3中，ψk(fi)根据重叠阵元之间的相位偏差计算得到，表示第k、k+1快拍之间由介质和路径扰动引起的频率fi的相位扰动。这正是ETAM被动合成孔径的关键所在对不同频率的相位扰动进行实时补偿，以实现孔径合成。ρm(fi)表示第m阵元频率fi的补偿因子，按照公式(4)计算得到。
实施例 下面结合某次海试数据和附图对本发明的具体实施方式
做进一步的详细描述。
试验参数拖曳阵水听器数目M＝40，水听器间距d＝1m；拖曳阵做匀速直线运动，速度v＝3.2m/s，信号采样率fs＝2000Hz。目标频带范围130～190Hz，方位125度，声速c＝1516m/s，快拍长度N＝2048。
需要注意的是取的数据长度越长，数据矩阵维数越大，运算速度会相应减慢。为了保证运算速度能够满足实时处理的要求，数据快拍长度不宜过大。一般FFT运算点数在2048点及以下均能满足要求。
采用基于ETAM被动合成孔径技术的宽带频域快速波束形成方法，具体步骤如下 步骤1对应图3中的301、302、303，用拖曳阵接收空间信号，取第k、k+1快拍的信号在时间域上做2048点FFT运算，得到快拍信号在不同频率的响应。如下式，N表示快拍长度，行表示时间采样，列表示阵元。
这里两个快拍之间的时间间隔τ取为拖曳阵匀速直线运动

个阵元距离所用的时间，依照“发明内容”中公式(1)可得 如图1-1所示，这样保证第k快拍、第21～40阵元和第k+1快拍、第1～20阵元在空间上重叠，为计算相位扰动做准备。
步骤2对应图3中的304，根据k、k+1快拍之间重叠的20个阵元，计算目标频带范围内的相位扰动。
需要说明的是，目标频带范围是130～190Hz，采样率2000Hz，FFT长度2048点，那么目标信号对应的离散频点是(130～190)/2000*2048＝(133～195)点，因此整个过程只需要对这一段频率范围进行处理即可。
相位扰动计算公式如下 m表示阵元序号，fn表示离散频点，按照前面的分析结果，n的取值范围是133～195。ρm(n)表示第m阵元频点n的补偿因子 这里Q是补偿因子的频带宽度，取为10。
步骤3对应图3中的305，对快拍k+1，第21～40的阵元数据进行相位补偿，实现孔径合成。
具体做法是对每个阵元不同频点的响应Xk+1，m(n)分别作ψk+1(n)、ψk(n)、Λ、ψ2(n)的相位补偿，如下式所示，将第k+1快拍、第20+m阵元虚拟到第1快拍、第40+20*(k-1)+m阵元。
步骤4重复步骤1、2、3，进行K次孔径合成，得到40+20*K个阵元的频域信号XETAM，fmin～fmax表示目标信号对应的频点范围。
步骤5对应图3中的306，对ETAM被动合成孔径得到的阵列频域信号在空间域上补零，延长阵元长度至MP，再对其在空间域上做FFT运算；将获得数据的零频分量移至谱中心。

式中，行fmin～fmax表示目标信号的频点范围，列1～MP表示某频点信号在空间域做FFT运算后的序号。
空间域上补零的目的是保证阵元数目远大于波束数目，使得空间域FFT运算结果和阵列波束输出的对应关系更加均匀，波束输出也就越精确。空间域补零可以有多种选择。既可以直接在原数据后直接补零，也可以在原数据中插值补零，只要补零后的数据长度满足要求即可。本实例直接在原数据后延长补零。
将零频分量移至谱中心是一种常用的处理方法，在许多文献中都有说明，对本领域技术人员来说，理解和实现该处理是可以胜任的。
步骤6对应图3中的307、308，在每一个频点上依据频率——波数网格对空域FFT运算结果进行校正，得到不同频点的波束输出矩阵L，其行表示搜索方位，列表示不同的信号频点。

步骤7对应图3中的309，将步骤6获得数据矩阵L的每一行的共轭转置与其自身相乘，结果即为在搜索方位上的输出功率。
P＝L*LT＝[P(θ1) P(θ2) Λ P(θ180)]T 步骤9重复以上步骤，进行下一次被动孔径合成波束形成，得到下一个时刻的目标方位功率输出。
图4对比了采用本发明前后，海试数据的波束输出。红色点划线表示采用本发明得到的波束输出，合成孔径后阵元数目200；蓝色虚线表示直接对40阵元实孔径信号做频域宽带波束形成的结果。从图中可以清楚看出，采用本发明得到的目标波峰的半功率波束宽度是未采用本发明的三分之一，也就说本发明能够有效提高目标方位的分辨率。
图5对比了采用本发明和对同样长度实孔径阵列的波束处理结果。红色点划线表示采用本发明，根据20阵元合成40阵元后的波束输出；蓝色虚线表示直接对40阵元实孔径信号做频域宽带波束形成的结果。对比可以看出，两种处理效果基本一致，也就是说本发明完全可以利用小孔径阵列得到和实际大孔径阵列同样的方位分辨率和作用距离。需要说明的是，获得这种好处存在一定前提在一次合成孔径时间内，阵列必须保持匀速直线运动。
为了说明本发明的实时处理优势，专门统计了matlab程序处理海试数据需要的时间长度。信号采样率fs＝2000Hz，数据快拍长度为2048点，直接对40阵元实孔径信号做频域宽带波束形成大概耗时0.06s；采用本发明大概耗时0.21s，合成孔径后阵元数目200。实际上，长度为2048点的快拍数据对应时间长度为2048/fs＝1.024s，可以看出，采用本发明完全可以实时处理海试数据。
图6是直接对40阵元实孔径信号做频域宽带波束形成的方位历程图，图7是采用本发明得到的海试数据方位历程图，可以看出，本发明得到的目标波束非常精细，更加有利于检测远距离微弱目标。
总之，本发明能够快速实时得到高分辨的目标方位功率输出。
最后所应说明的是，以上实施例仅用以说明本发明的技术方案而非限制。尽管参照实施例对本发明进行了详细说明，本领域的普通技术人员应当理解，对本发明的技术方案进行修改或者等同替换，都不脱离本发明技术方案的精神和范围，其均应涵盖在本发明的权利要求范围当中。
权利要求
1.一种快速宽带频域扩展拖曳阵波束形成方法，该方法将ETAM被动合成孔径技术和宽带频域波束形成方法结合在一起，可用于对微弱目标的检测和定位，同时速度较快，便于实时处理；所述的方法包含如下步骤
1)用线阵接收空间信号，得到M个阵元的时域信号；取重叠阵元数目为N，其中N大于零且小于阵元数目，计算两次快拍数据的时间间隔阵τ，即阵运动M-N个阵元距离需要的时间；
2)对第k快拍和第k+1快拍数据在时间域上做快速傅里叶变换；
3)确定目标辐射信号的频带范围[fmin，fmax]，根据第k快拍和第k+1快拍的重叠阵元计算相移参数ψk+1(fi)
k表示快拍序号，n表示阵元序号，fi表示频率，ρn(fi)表示第n阵元频率fi的补偿因子
4)对第k+1快拍、第M-N+1～M阵元信号的不同频率分量分别进行相位补偿，将其虚拟成第k快拍的第M+1～M+N阵元；依次类推，这N个阵元信号可以虚拟成为第1快拍的第M+N*(k-1)+1～M+N*K阵元；
5)重复以上步骤，进行K次的孔径合成，可以虚拟得到MM＝M+N*K个阵元接收的信号在[fmin，fmax]范围内的频域响应；
6)对MM个阵元信号在空间域上补零，再对其在空间域上做FFT变换，然后将获得数据的零频分量移至谱中心；
7)在处理带宽内的每一个频点上依据频率——波束网格对其进行校正，得到不同频点的波束输出；
8)合成不同频率的波束形成结果，得到MM个阵元在不同方位上的功率输出，从而实现目标的精确定位。
2.根据权利要求1所述的快速宽带频域扩展拖曳阵波束形成方法，其特征在于，步骤1)所述的重叠阵元占阵元总数的一半。
3.根据权利要求1所述的快速宽带频域扩展拖曳阵波束形成方法，其特征在于，所述步骤6)中，采用FFT进行波束形成，所述的波束输出公式(7)如下；
所述的阵列信号的空间域FFT运算公式如下
所述的波束方位θs和FFT变换域k关系如下
其中，k是FFT的运算序号，fi表示频率，d是阵元间距，θs是扫描方位取值范围0～180度，MM是经过ETAM合成孔径运算后的阵元数目。
4.根据权利要求1所述的快速宽带频域扩展拖曳阵波束形成方法，其特征在于，所述的空间域FFT运算结果用频率——波束网格对其进行校正得到单频信号的波束输出；
所述的频率——波束网格校正是预先产生的，计算公式如下所示，floor表示取整运算；
5.根据权利要求1所述的快速宽带频域扩展拖曳阵波束形成方法，其特征在于，步骤6)所述的空间域补零，具体采用直接在原数据后直接补零，或在原数据中插值补零。
全文摘要
本发明提出一种快速宽带频域扩展拖曳阵波束形成方法，该方法将ETAM被动合成孔径技术和宽带频域波束形成方法结合在一起，包含如下步骤1)用线阵接收空间信号，得到M个阵元的时域信号；取重叠阵元数目为N；2)对第k快拍和第k+1快拍数据在时间域上做快速傅里叶变换；3)确定目标辐射信号的频带范围，根据相邻任意两快拍的重叠阵元计算相移参数4)对不同频率分量分别进行相位补偿得到虚拟的合成阵元；5)重复以上步骤，进行K次的孔径合成；6)对合成的阵元信号在空间域上补零做FFT变换，然后将获得数据的零频分量移至谱中心；7)在每一个频点进行校正，得到不同频点的波束输出；8)合成不同频率的波束形成结果，从而实现目标的精确定位。
文档编号G01S7/539GK101813772SQ20101016326
公开日2010年8月25日 申请日期2010年4月30日 优先权日2009年12月31日
发明者李峥, 薛山花, 田彪, 李宇, 黄海宁 申请人:中国科学院声学研究所, 北京中科海跃科技有限责任公司