基于被动声信号解算方法的高速目标轨迹测量

2023-04-03孙贵新于文彬

计算机测量与控制 2023年3期

孙贵新,于文彬

(1.中国人民解放军91550部队230所，辽宁大连 116023;2.上海交通大学电子信息与电气工程学院，上海 200240)

0 引言

脱靶量测量是武器系统鉴定中不可缺少的一项重要指标，是外场试验的核心测试参数[1]。尤其对于新一代外场，高速目标呈现出命中精度不断提高、速度不断提升、射速不断加快、全天候等趋势，这给高速目标脱靶量的检测评估带来极大挑战，主要体现在3个主要方面：需要更高的目标脱靶量检测效率、需要更准确全面地测量目标攻靶信息、需要满足全天候测试的要求。

常用的脱靶量测量技术方法有：靶板法、无线电技术[2]、天幕靶技术[3]、光学测量方法[4]、GPS测量方法[5]等技术方法。靶板法即传统的人工测量方法，存在着准确性低、强度大、实时性差、效率低、危险性高且无法测量连续多个目标的不足；虽然光学、电磁波都能对目标攻靶过程的末段飞行和落点进行观测，每种技术手段各有优势和短板，没有任何一种单一手段能够实现对高速目标攻靶的全天候高精度观测评估，具体来说，基于光学的测量方法在一定程度上，系统具有组成相对简单、成本较低、测量精度较高、分辨率高等优点，但很难适用于夜间和雨雪天气，而且无法满足大范围、远距离场景下对精度的要求；对于无线电技术，虽然其定位精度高，自动化程度高，但只适用于特定范围，例如由于地面或者海面反射波的干涉作用造成盲区，从而使无线电技术失效；红外和微波等观测手段受雨雪等天气的影响也很严重。而声学场方面，由于对电磁干扰、雨雪及夜间能见度等因素不敏感，且观测范围大、成本低廉、体积较小等优点[6]。因此，可满足在不同环境下对不同种类和目标的动态飞行轨迹脱靶量的高精度和高准确率的测量观测，弥补恶劣环境下观测维度不足的劣势。

基于声学原理的声探测技术，是利用声学基阵接收目标的声特征，利用方位、时延等信息实现目标的脱靶量测量，因其特有的优点，现阶段存在大量基于声学原理的脱靶量测量研究成果[7]。文献[8-9]中根据超音速目标特有的激波到达时间，提出基于时间梯度和波前方向矢量的目标测量方法，实现超音速目标位置估计的同时，给出了速度的估计方法，最后通过仿真实验验证算法的有效性，但该算法的有效性是在较大范围内布设多个测量基站(基站间间隔为80 m)为前提，这在实际较小几何尺度的靶船试验场景难以实现；文献[10-11]中利用声学测量受天气等因素影响小的优势，利用声学基阵采集的声特征实现目标测量，弥补预警雷达失效或精度低时的目标测量，但该方法存在解算模型不明确以及解算精度未进行理论分析的不足[10-11]；另一方面由于声学基阵是保证声学测量精度的基础，因此，对声基阵的研究十分必要，对此研究人员根据不同测量需求及外场场景设计了四元阵[12-13]、五元平面阵[14]、五元空间阵[15]、正四棱锥阵[16]、七元空间声阵[17]、T型靶[18]等阵型，并在此基础上给出相应的解析解，但在实际外场试验中，场景或者待测量目标往往并不唯一，因此不具备广泛适用的能力。

针对现有方法的不足，结合高速目标和新一代外场试验场景特殊要求，本文提出了基于声学的轨迹测量方法，实现高速目标轨迹跟踪，仿真试验验证了对超音速目标其测量精度优于10 m，对高亚音速目标测量精度优于10 m。提出的测量框架弥补了外场观测手段的空缺，同时极大程度丰富了外场观测维度。

1 系统架构

声学轨迹测量系统主要由轨迹测量子系统、数据采集与存储子系统、地面显控子系统3个部分组成，实现了各类目标信号高采样率采集、处理和存储，前端探测节点状态等信息实时更新，事后轨迹解算结果的直观显示。系统框架具体如图1所示。

图1 声学轨迹测量系统架构图

如图1所示，对于轨迹测量子系统，其中各声学基阵用于采集目标经过空气产生的摩擦声，采集的声信号也会因目标的外形、速度等不同而具有不同的特征，采集的声信号用以作为目标轨迹解算的数据源；气象测量单元用于测量不同任务场景下的气象数据，如风向、风速、温度和湿度等，实现对环境的准确描述以提升目标测量精度；信号采集与处理单元，通过信号调理模块实现多路声信号的滤波、放大、AD(analog to digital, AD)转换，其中滤波是通过Butterworth带通滤波器完成带外噪声的滤除以及外界电磁高频噪声，AD模块实现将采集的声模拟信号转换为电信号，其转换精度影响目标测量精度，综上，通过信号调理模块完成信号的采集与转换，去除或分离声信号中所含背景噪声或干扰分量。数据采集与存储子系统对各探测节点预处理数据、状态信息、测量信息的读取、存储和转发，其中对于存储介质的选取，需要具有存储容量较大、尺寸较小、防震防摔等优点，以提高系统在工作时的环境适应能力。地面显控子系统内含声学轨迹的解算算法和显示界面，实现地理信息的三维可视化显示，目标轨迹的实时显示，节点状态与相关信息的查看，而且，地面显控子系统包含通信基站，实现相关指令的下达和数据抓取，另外，地面显控子系统实时监测系统状态及异常行为，试验或仿真测试中所有操作及命令记录日志，为事后查看相关实验流程提供方便。

2 实现过程

2.1 目标轨迹测量声学基阵

目标轨迹测量声学基阵主要由一定数量的声学基阵按照一定空间几何形状排列而成，主要参数包括各声学基阵内部参数、声学基阵数目、声学基阵间间距。目标轨迹测量声学基阵的阵型设计原则是，根据外场试验场景及测量任务，采用最少的声学基阵实现较高的目标测量精度。本文拟采用五条声学基阵，具体排布结构示意图如图2所示。声学基阵S0、S1、S2、S3和S4五条阵按图中几何位置布放，对应的直角坐标分别为S0(0,0,0)、S1(D,0,0)、S2(0,D,0)、S3(-D,0,0)和S0(0,-D,0)，其中D为声学基阵阵间的距离。

图2 声阵列布阵示意图

对于每个声学基阵，是由声传感器阵列组成，声传感器阵列主要参数为声传感器数目、阵列孔径大小、声传感器阵元间距及其空间分布形式等。根据外场场景尺度不同和测量任务不同，常见的声传感器阵列结构主要有线列阵、平面阵和立体阵等。其中，线列阵结构简单，设计容易，计算难度小，但目标的空间定位效果不佳，且在没有先验信息的情况下具有左右舷模糊的不足。平面阵列计算量适中，可用于对目标的平面和空间定位；立体阵列空间定位效果良好，但结构较为复杂，且计算量较大。在本文中，考虑到外场试验场景、任务需求，以及硬件系统主要基于FPGA (field programmable gate array, FPGA)完成，且方位角和距离因素对目标测量的精度影响更为重要，故从硬件和环境因素出发，选用平面五元阵，可实现目标的平面和空间定位。

声传感器阵列间距D是由目标声信号的线谱频率、波长关系等确定，计算公式为:

(1)

式中,λ为声波波长，其取值范围为0.085～3.4 m；c为声音传播速度；f为信号频率。

若假设声信号在空气中的传播速度为340 m/s，为确保不出现目标模糊，可得声传感器阵元间距的取值应小于1.7 m。可根据目标声信号的主频率与基阵阵列间距D的关系式确定阵元间距D。总体上说，阵列的孔径越大，定位效果越好。具体应用时，需综合定位功能、便捷性和性价比等综合考虑。

声阵列装置包括支架、一个或多个直线传声器阵列单元及其附件等。其中，传声器阵列单元包括圆管及其安装座、传声器安装座、传声器和防风球等，属于结构可调伸缩式设计模式。具体来说，传声器固定在传声器安装顶部中心孔内，传声器安装座底部与圆筒顶部同轴螺纹相连，圆筒底部通过圆筒安装座进行安装。同时，支架外表面包裹吸声材料，且高度可调。可综合任务场景地质条件、气象条件、运输、及可操作性等因素进行选材和结构设计。同时，为确保声阵列长时间户外工作的可靠性，还将在声阵列上配备安装避雷针和浪涌保护器，以保护目标轨迹测量系统在户外能长时间的正常运行。轨迹测量声阵列示意图如图3所示。

图3 轨迹声阵列示意图

2.2 信号采集与处理

信号采集与处理是本系统中用于完成超音速和高亚音速目标的脱靶量测量的第一阶段的工作，即目标在空中飞行时与空气摩擦时产生的声信号的采集和预处理，其主要包括顶层软件模块、声信号采集软件模块、声信号预处理模块、数据存储模块等，信号采集与处理结构示意图如图4所示。

图4 信号采集与处理示意图

由图4可知，信号采集与处理中的顶层模块是整个声学轨迹测量系统的入口；声信号采集软件模块通过动态信号采集卡完成对声信号的采集，然后经工业硬盘完成对采集数据的存储，供后续分析处理使用；声信号预处理模块基于内嵌算法对采集的声信号中的干扰分量及背景噪声等进行滤除或分离，下面进行具体介绍。

2.2.1 信号采集软件

信号采集工作主要在采集卡上实现，工作状态下的数据采集软件流程图如图5所示，此处，首先动态信号采集卡处于等待状态，当接收遥控指令后，开始工作完成初始化及状态自检，无误后等待接收外部触发信号，接收到触发信号后进行声信号的存储及预处理。

图5 数据采集流程图

2.2.2 声信号处理

对于目标在最后攻靶段与空气摩擦产生的声信号，声学基阵采集到的是带噪声信号，为实现背景噪声或干扰分量的去除或分离，本文拟采用滤波、经验模态分解(empirical mode decomposition, EMD)算法和小波阈值相结合的处理方法，其中滤波是由硬件系统完成，即通过Butterworth带通滤波器完成指定频带外的噪声滤除以及外界电磁高频噪声。

对于EMD算法，其核心是将非线性、非平稳的原始信号分解为一系列表征原信号某一频带或某一尺度的本征模态函数(intrinsic mode function, IMF)的组合，令瞬时频率在IMF上存在一定的物理意义。另一方面，由于EMD分解具备正交性和完备性的数学性质，使得分解得到的一系列IMF分量不仅能够保留原始信号的特征信息，并可以按一定准则重构原始信号。因此，相比于直接处理原始的非线性、非平稳信号，处理经分解得到的一系列IMF分量，在提高计算效率的基础上，还能够避免直接处理原始信号造成的信息缺失及特征分析中的人为因素干扰，提升原始信号去噪的客观性。

2.2.3 数据采集与存储

靶区数据采集与存储子系统，用于采集和存储轨迹测量声阵列、目标轨迹测量声阵列和气象测量单元的状态信息、量测信息等；同时将各个探测节点的状态信息(含全球导航卫星系统(global navigation satellite system ,GNSS)时间位置信息)等的存储和转发，因此据采集与存储子系统具备声信号、震动、气象等测量数据及各个声学基阵的探测节点预处理数据、状态信息、测量信息的读取、存储和转发功能。

对于各探测节点测量数据、状态信息和测量信息，根据“容器化+微服务”架构要求，进行数据清洗和整合，将各非结构化的原始数据转化为结构化数据，方便后续数据的快速准确调用。在存储数据前需要进行数据预处理，具体来说包括对震动数据、气象数据的预处理，以及目标声信号原始波形滤波、信号检测、直流偏置去除、插值与拟合等。

2.2.4 地面显控

地面显控主要是根据处理后的声信号实现目标轨迹解算，并实现解算结果的实时显示。超声速或高亚音速目标的轨迹声信号解算包括两个部分：单个声学基阵接收的声信号对目标声源的时变方位信息解算或声强极值点解算、和基于空间信息的目标运动轨迹解算。超音速目标和高亚音速目标轨迹解算示意图分别如图6和图7所示。此处，以3个声阵列为例进行说明。

图6 轨迹声信号解算流流程图(超音速)

图7 轨迹声信号解算流流程图(高亚音速)

对超音速飞行的目标来说，其处理思路为利用各个声阵列接收的信号实时解算目标相对于阵列的方位，然后根据不同阵元获得的方位量测值，进而利用多面交汇的方式获得目标运动轨迹。

对高亚音速飞行的目标来说，主要是利用各个声阵列接收的声强极值点来分别估计目标运动轨迹的垂线方向，通过计算多个垂线的向量积实现对目标运动方向的估计，进而利用多面交汇的方式获得目标运动轨迹。

对于方位估计，根据测向原理，传声器阵列的测向算法可分为时延声源测向、波束形成声源测向和超分辨空间谱估计声源测向。其中，超分辨空间谱估计声源测向能够突破瑞利限制，且对阵列构造要求不太严格。多重信号分类(multiple signal classification, MUSIC)算法由Schmidt等提出，通过数学矩阵变换将信号分为信号子空间和噪声子空间，然后利用这两个空间的正交性来估计信号的参数。该方法通过在全角度空间进行谱峰搜索来估计声波方向，有着较高的测向精度，抗噪能力强，且可实时处理和能对多个信号进行测向等。本文拟采用2D-MUSIC算法来实现目标声源方位信息的估计。

对于阵列间距D的M个阵元来说，若入射信号源为S，方位角为θ，俯仰角为φ，且假定信号子空间和噪声子空间正交时，MUSIC算法的谱估计公式为:

(2)

式中,N为声信号的频域快拍数；UN为声信号周围环境噪声的噪声子空间。

改变θ和φ变化，通过搜索谱峰来估计目标声源的方位角和俯仰角。

在实际应用中，由于接收到的声信号具有多个频率分量，需要首先将声信号在频域上分解为多个分量。将所有分量经过式(1)得到的结果进行加权，即可得到声源信号最终的空间谱。具体过程如下：

对每个分量计算其对应的谱密度矩阵，为:

(3)

计算第j个子带对应的二维空间谱函数，为:

(4)

进行加权平均，构造声源信号的二维MUSIC空间谱函数，为:

(5)

在实际应用中，由于接收到的声音信号具有多个频率分量，需要首先将声音信号在频域上分解为多个分量。将所有分量经过式(1)得到的结果进行加权，即可得到声源信号最终的空间谱。

理想情况下，即每个声阵列测得的方位角和俯仰角都不存在误差，则3个声阵列的方向线在目标声源处会相较于一点，此时只需要求解三条方向线的交点就可以获得目标声源位置。

但在实际应用中，由于侧向误差的存在，三条方向线在空间中可能互为异面直线，不存在直接交点，如图8所示。图中，S为目标声源；L1、L2和L3为三条方向线；S′为目标声源位置近似估计。

图8 空间异面直线等效交叉示意图

3 算法介绍

3.1 超音速目标轨迹解算

在一定范围内按一定几何形状布设多个声学传感器形成声学基阵，各个声学传感器在接收超声速目标的激波信号后，求解各时刻超音速目标相对声学基阵的方位，最后根据平面相交原理实现对超音速目标的矢量脱靶量测量。

马赫波是指当目标超音速运动时，对运动周围的空气形成扰动。马赫波叠加并转变为激波，如图9所示。

图9 目标激波示意图

图9中的θ为马赫角，其数学关系可描述为：

(6)

其中:Ma是马赫数。

当目标运动速度大于声速时，满足Ma>1，从而构成马赫锥，当目标运动速度越大时，马赫角越小。

另一方面，激波的生成是指在目标运动时，目标的尾部产生真空区，其示意图如图10所示。

图10 激波压力曲线

其中:τ1、τ2及T是时间间隔。

如图10所示，超音速目标运动时首先产生马赫锥，而马赫锥进一步产生的凹角，令压力从常压P0迅速升到P0+P1，紧接着马赫锥产生的凸角，令压力从P0+P1迅速降到P0，在此过程中的柱面膨胀，令压力从P0迅速降到P0-P1，超音速目标底部变化令压力从P0-P2快速升到P0。

图11 超声速目标过靶示意图

图12 目标运动轨迹与声基阵的几何关系示意图

(l1,m1,n1)·(l,m,n)T=cos(90°-θ)

(7)

同理5个声基阵可构建5个观测方程为:

(8)

通过式(8)可知为齐次方程，其中未知数为(l,m,n,θ)。

在试验中，目标的飞行速度可通过遥测数据获得，通常为已知量，故可利用马赫角求解出θ，因此，式(8)可改写为:

(9)

式(9)中，(li,mi,ni)为观测值，sinθ为已知值，可求解出未知数(l,m,n)。

3.2 高亚音速目标轨迹解算

对高亚音速运动的目标来说，声传感器接收到的声信号持续时间很短。在如此短的接收声信号中，难以将其中某一段信号对应到相应的目标轨迹点上，故在确定该类飞行目标的运动轨迹时，进行如下基本假设：

1)声强最大点是目标运动轨迹与声阵列中声传感器最近的点；

2)目标运动轨迹直线位于基准传声器接收的声强最大点入射波矢量垂直的平面上；

3)目标运动轨迹与矢量所在的平面垂直。

根据基本假设，利用两个以上的入射波矢量，即两个以上的超短基线阵，可估计出目标运动轨迹的方向向量，结果可表示为:

(10)

目标运动方向向量归一化估计值为:

(11)

(12)

利用多个平面方程相交，可以得到目标运动轨迹上某一点的坐标，从而完成对目标运动轨迹的估计。

4 仿真分析

首先，给出超音速目标和高亚音速目标轨迹测量方法及步骤。

4.1 超音速目标轨迹解算

首先根据外场试验场景和外场试验环境，合理假设目标的运动形式和运动参数，声学基阵的位置和声学基阵间间距，合理假设声学基阵和目标位置和相关参数的误差，在此基础上可获得仿真的声学基阵测得的关于方位的量测数据；各声学基阵获得的声信号根据式(1)～(4)的原理解算目标的方位，然后通过平面相交原理及其误差情况下无法交汇于一点的补偿方法(示意图8)，获得目标的位置，最后根据目标真实位置进行算法的目标轨迹测量精度分析。

4.2 高亚音速目标轨迹解算

对于高亚音速目标的轨迹解算，其思路与超音速目标轨迹结算思路相同，不同点在于高亚音速目标在仿真中，量测数据为声强，即根据式(10)求得目标轨迹的方向向量，然后通过最小二乘估计求得平均值，在此基础上，通过平面相交原理求得目标轨迹，最后根据目标真实位置进行算法的目标轨迹测量精度分析。

下面进行具体介绍，外场声学基阵仿真示意图如图13所示。

图13 外场声学基阵仿真示意图

根据图11设置靶船外场试验仿真场景图13，在图13中，构建载体坐标系，其中声学基阵为S1-S5按图中所示几何结构布放，其中基阵间间距根据靶船几何尺度进行合理假设，对于目标运动轨迹，其运动形式假设为匀速直线运动，当目标经过靶面时，各声学基阵采集信号，然后通过处理，实现目标轨迹解算。

具体来说，根据上述外场试验示意图设定仿真参数，即5个声学基阵的中心坐标分别设为(-35,-15,2)、(-35,15,2)、(0,0,2)、(35,-15,2)、(35,15,2)(单位：m)。假定目标航向角为90°，俯仰角为3°，速度v=2.5马赫，各个声学基阵测向角误差精度为0.1°。在过靶区域每隔1 m设置一个过靶点，分别计算在不同过靶点上的误差，以过靶点坐标为(0,0,15)为例，过靶点定位结果为(0,0,14.89)，当武器穿过靶面时，靶船首尾方向误差10 m，高度方向定位结果为14.89 m，定位误差为0.11 m；武器航向角误差0.11°；武器俯仰角误差0.28°。过靶点位置、z轴方向、x轴方向、俯仰角、航向角仿真结果如图14所示。

图14 超音速目标脱靶量结算误差仿真结果

由图14可知，当目标与靶船航向垂直时，脱靶量解算精度最高，平行时，解算精度较差；当目标水平攻击靶船时，脱靶量解算精度最高，垂直攻击时，解算精度较差；目标过靶点靠近靶船首尾时高度方向误差呈增大趋势；随着目标过靶点高度降低,靶船首尾方向误差呈增大趋势。算法的过靶点解算总体精度优于10 m。

声学基阵的中心坐标分别设为(-35,-15,2)、(-35,15,2)、(0,0,2)、(35,-15,2)、(35,15,2)。过靶区域-35≤x≤35，5≤z≤30，y=0。假定目标航向角为90°，俯仰角为3°，速度v=0.75马赫，声学基阵各向测角误差精度为0.1°。在过靶区域每隔1 m设置一个过靶点，分别计算在不同过靶点上的误差。过靶点位置、z轴方向、x轴方向、俯仰角、航向角仿真结果如图15所示。