周家新, 陈建勇, 单志超, 陈长康
(1.海军海洋测绘研究所, 天津 300061; 2.海军航空工程学院 电子信息工程系, 山东 烟台 264001)
铁磁性物体受到地磁场磁化产生的磁场,与环境原有的背景磁场叠加在一起,使得铁磁性目标所在区域的空间磁场分布发生改变,从而形成磁异常。通过识别这种磁异常信号进行目标检测与定位的技术称为磁异常检测[1-2]。由于航空磁探仪受空气、海水、泥沙土壤等介质的影响小,且具有可以连续搜索、使用简单可靠、定位精度高、反应迅速、搜索效率高等优点[3],作为被动检测技术的一种,已成为水下潜艇目标的重要探测手段,被广泛应用于地磁勘探、管线搜索以及未爆弹药探测等[4-6]领域。
目前,磁异常检测技术主要分为两类:一是基于目标特征,将铁磁性目标视为磁偶极子模型,从而给出目标磁异常的数学表达式,进而使用统计分析理论进行目标信号检测,主要有正交基函数(OBF)磁异常检测器与主成分分析(PCA)检测器[7-10];二是基于磁背景噪声的统计特性,通过噪声与目标信号之间的差异,借助统计学中的假设检验进行判决,主要有最小熵检测器(MED)与高阶过零(HOC)检测器[11-14]。
传统基于目标特征的OBF磁异常检测器[7]和PCA检测器[10]通常基于Anderson函数将目标信号进行正交基分解,得到能量信号进行阈值检测。该方法需要对检测平台磁补偿后剩余磁噪声的特性进行分析和处理,使其满足能量正交分解最优检测的相关要求。基于噪声特性的MED磁异常检测通常将磁背景噪声视为高斯分布,但在复杂磁环境中,磁噪声分布具有非高斯性,主要的磁噪声源为地磁噪声[13],其功率谱密度为1/fα,其中f为地磁噪声频率,指数α的范围为0~2. OBF检测算法被证明是高斯白噪声情况下的一个最优检测方法,但该方法依赖于背景噪声特性。而在实际磁场测量环境中,噪声功率谱密度为1/fα更常见。另外,MED检测器[11]与HOC检测器[12]需要大量的先验噪声样本,在低信噪比的实际应用中往往难以满足检测需求。
文献[15]提出了一种适用于非高斯噪声背景的带通滤波结合OBF分解的磁异常信号检测算法,将非高斯噪声进行白化处理,基于OBF分解得到能量信号并用于检测铁磁性目标。但由于航空磁探潜的背景噪声干扰强,该检测算法的性能有限。文献[16]提出了一种航空磁探中水下目标的自适应探测方法,基于学习飞行补偿背景磁场,并对噪声进行自适应白化处理,从而对目标信号进行检测。但该方法的学习飞行不稳定,对背景磁场补偿效果有限,难以将目标信号从强磁背景干扰中准确识别。
为提高航空磁探仪在低信噪比条件下对水下铁磁性潜艇目标的检测概率,降低航空磁探潜过程中的虚警概率,本文提出了一种联合估计检测方法。基于磁偶极子建立潜艇目标磁场信号模型,使用航空磁力仪观测样本,通过遗传算法(GA)在潜艇目标先验信息范围内搜索目标的磁矩矢量、位置矢量和航向等参数,初步检测目标是否存在并得到估计信号。使用目标的估计信号与测量平台的实测信号进行平均相关检测,从而更进一步确认潜艇目标的存在,得到目标的位置、磁矩等参数,同时实现目标的检测、定位和分类。
航空磁异常探潜中,航空磁探仪和水下潜艇目标的距离一般都在数百米之遥,甚至更远,根据文献[17]的研究结论,在这种情况下,目标磁场可以等效为偶极子磁场[18-19]。
对于一般的铁磁性潜艇目标,其在空间点P(x,y,z)处产生的磁场强度H可以根据磁偶极子模型计算得到,以水下潜艇目标中心为坐标原点(见图1),则可以得到P点处的磁场强度为
(1)
式中:M为潜艇目标的磁矩矢量;r0为目标源点O指向空间测量点P的单位矢量;r为目标源点O与空间测量点P之间的距离;o(1/r5)为高阶次极子在铁磁性潜艇目标磁场中的贡献。随着测量点和目标之间距离的增大,高阶次极子的作用逐渐变小。在忽略高阶次极子的作用之后,水下铁磁性潜艇目标磁场信号模型可表示为
(2)
图1中潜艇目标磁矩M可以表示为3个分量的形式,分别为Mx、My、Mz,将磁偶极子模型转换到直角坐标系下,磁偶极子模型在空间产生的磁场矢量可以表示为H=ixHx+iyHy+izHz,ix、iy、iz为3个坐标轴的单位矢量,则三分量Hx、Hy、Hz为
(3)
在航空磁异常探潜中,将磁场强度H在真空(空气)中引起的磁感应强度记为B,并有如下关系:
(4)
式中:μ0为真空磁导率,μ0=4π×10-7H/m. 在海水中磁导率近似等于真空磁导率,而且由于磁感应强度的基本单位T非常大,航空磁探测中多使用nT作为单位。
将(3)式代入(4)式,得到水下潜艇目标在空间产生的磁感应强度三分量Bx、By、Bz为
(5)
使用(5)式计算时,右边仍采用国际单位制,左边所得结果单位为nT.
将潜艇目标空间三分量磁场转换为航空光泵磁探仪测量得到的空间磁场总场信号,设θ为磁北方向到潜艇目标航向的逆时针角度(航向角),将目标的三分量磁感应强度投影到地磁场方向,则目标磁感应强度的水平分量Bh和垂直分量Bv的形式为
(6)
沿地磁场方向的单位磁场向量为(cosγ,sinγ),其中γ为测量点处的地磁倾角,则铁磁性潜艇目标的光泵测量值即为目标磁场向量沿地磁场方向的投影:
B=Bhcosγ+Bvsinγ.
(7)
将(5)式、(6)式代入(7)式中,整理得到铁磁性潜艇目标光泵测量信号模型为
(8)
目标磁感应强度三分量Bx、By、Bz及其投影到地磁场方向的光泵总场信号如图2所示,其中该信号以磁偶极子为中心,测线中心位于目标中心正上方,高度为200 m,测线的范围x∈[-450 m,450 m],磁矩为Mx=80 000 A·m2,My=-20 000 A·m2,Mz=40 000 A·m2,目标航向角为0°,地磁倾角为45°.
由(8)式可知,铁磁性目标信号有8个参数,若其中地磁倾角γ可以认为已知,则仍有位置x、y、z,磁矩Mx、My、Mz,航向角θ等7个未知参数。
航空磁探仪通过在任务区域的飞行来实施航空磁异常探潜任务,建立如图3所示的探测坐标系。
经典探测宽度模型[20]可以表示为
(9)
式中:Rm表示航空磁探仪的作用距离;Hp表示磁探仪所在高度;Hs表示目标所在深度。设航空磁探仪直线飞行,从磁探仪采样信号开始到结束的飞行距离为L,为简化探测过程,将探测区域近似为一个长方体,从而可以得到目标所在区域长为L+W,单位均为m. 基于磁探仪的探测宽度模型和磁偶极子信号模型,能够实现磁性目标的远程定位,有利于快速缩小目标的搜索范围。
航空磁探潜过程中,磁力仪采样一段信号,信号长度为N,观测信号分别记为B1,B2,…,Bn,…,BN,即
(10)
由(8)式的信号模型可知,待估计参量矢量为
a=[x,y,z,Mx,My,Mz,θ]T.
(11)
根据文献[21]设置探测区域以及待估计参量矢量的范围。基于目标先验信息,在探测区域内,使用GA对待估计参量进行搜索[22-24],基于Fréchet距离构造适应度函数。基于GA的参量估计使得适应度函数在搜索过程中总能收敛到一个极小值点,能以很大概率找到全局最优解。
Fréchet距离应用广泛[25],可描述两曲线之间的相似性程度[26]。在航空磁探潜水下目标的联合估计检测中,可以用于量化估计信号与实测信号的相似性。铁磁性目标的估计信号为B(),实测信号B(a)与B()为两个离散的点序列串,其中B(a)为〈B1(a),B2(a),…,BN(a)〉,估计信号B()为〈B1(),B2(),…,BM()〉,M为估计信号的长度。则Fréchet距离为
Fd(B(a),B())=max (Ed(B(a),B()),
min (Fd(〈B1(a),B2(a),…,BN-1(a)〉,
〈B1(),B2(),…,BM()〉),∀N≠1,
Fd(〈B1(a),B2(a),…,BN(a)〉,
〈B1(),B2(),…,BM-1()〉),∀M≠1,
Fd(〈B1(a),B2(a),…,BN-1(a)〉,〈B1(),
B2(),…,BM-1()〉),∀M≠1,N≠1)),
(12)
式中:Ed(B(a),B())为B(a)与B()之间的欧氏距离。根据文献[21]中的数据设定一组目标位置、磁矩及航向等参数,生成目标空间磁场分布,模拟磁力仪采样目标信号,得到仿真的实测信号B(a). 然后基于GA在区域范围内搜索待估计参量,搜索过程如图4所示,在规定的遗传代数(3 000)内搜索,若适应度未达到门限值则判定目标不存在。针对搜索过程,在空间直角坐标系中标记出估计参数和设定参数位置。图5为基于GA的潜艇目标在航空磁探仪搜索区域范围内的待估计参数搜索定位图,目标位置参数为x=500 m,y=100 m,z=150 m,磁矩参数为Mx=80 000 A·m2,My=-20 000 A·m2,Mz=40 000 A·m2,航向角参数为θ=190°. 在图5(a)和图5(b)中,蓝色方块为目标的真实位置以及真实磁矩,红色十字点为每次搜索定位的结果;图5(c)中,蓝色直线为θ=190°. 由搜索定位结果可以得到,基于GA的搜索目标参量,能够快速收敛到实际真值。
航空磁探潜中的磁异常检测属于双择检测,水下潜艇目标存在状态对应假设H0和假设H1两种情况。根据磁探仪采样的观测信号,水下潜艇目标双择检测问题表示如下:
(13)
式中:H0和H1对应两种假设情况,H0表示目标不存在,H1表示目标存在;x(n)为观测信号;ω(n)为噪声信号;n=1,2,…,N.
通过估计目标的磁矩矢量、位置矢量以及航向等参数,得到目标磁感应强度三分量的估计,进一步得到目标的估计信号。基于目标的估计信号与测量平台的观测信号进行平均相关检测,从而判定目标是否存在。构造如下检验统计量:
(14)
式中:输入x为观测信号矢量;RxB(m)为观测信号与估计信号的互相关函数。根据检验统计量构造的联合估计检测器实验框图如图6所示。根据(8)式的目标光泵测量信号模型生成目标的估计信号,通过相关器对估计信号与观测信号做平均相关,将得到的检验统计量与门限β进行比较,判定目标是否存在。检测门限β根据奈曼- 皮尔逊准则得到,由大量的观测噪声样本进行统计分析,计算得到最终的门限值。
采用光泵磁力仪实测的磁噪声数据和磁偶极子仿真目标信号数据进行联合估计检测器的性能分析,实测噪声数据片段如图7(a)所示。航空磁探仪飞行高度为150 m,飞行速度为100 m/s,飞行航向为0°,测量地区的背景磁场为49 699 nT,地磁倾角为52°,其中磁力仪的分辨率为0.01 nT,采样频率为100 Hz. 潜艇目标信号使用(8)式的磁偶极子模型仿真产生,目标信号加实测磁噪声经过去趋势处理后如图7(b)所示,目标出现在采样点2 500~2 800段。目标信号宽度为300,检测器输出中目标信号的宽度接近300. 为与文献[7]中的OBF检测器性能进行比较,将目标信号与实测磁噪声之间的幅度信噪比设置为0.3. 基于目标估计信号的联合估计检测器(ES-COD)归一化输出如图7(c)所示,其输出的幅度信噪比为10,能够大大提高检测器的输出信噪比。文献[7]中OBF检测器的输出幅度信噪比为5,相比之下,ES-COD检测器的输出幅度信噪比提升了1倍,检测性能有所提高。文献[15]提出的带通滤波结合OBF检测器分解的磁异常信号检测算法适用于非高斯噪声环境,但其输出幅度信噪比仅为2. 在压制噪声能力方面,与该算法相比,本文提出的ES-COD检测器性能更优,能够显著提高信噪比。
选取实际铁磁性潜艇模型代替潜艇目标对检测器性能进行Monte Carlo仿真分析,艇模三视图如图8所示,噪声信号采用实测噪声数据。实验中,在潜艇模型中心上方10 m的高度面上沿x轴方向每0.1 m取1个测量点,共计151个考核点,x∈[-7.5 m, 7.5 m]。在考核点处使用G858光泵磁力仪进行测量,磁力仪性能参数如表1所示。横向测量范围为15 m,潜艇模型长度为2 m,测量范围达到15倍艇长。测量高度为10 m,是艇长的5倍。测量结果如图9所示。基于奈曼- 皮尔逊准则,在虚警概率为0.001的条件下,各信噪比情况仿真10 000次,得到如图10所示的检测概率曲线。
铯光泵参数铯光泵数值测量范围/nT18000~95000分辨率/nT0.05温度漂移/(nT·℃-1)<0.05梯度容限/(nT·in-1) 500测量精度/%<1频域噪声/(pT·Hz-12)<4
分别将文献[3]、文献[7]和文献[12]中的最小均方估计- 正交基函数(LMS-OBF)检测器、OBF检测器和HOC检测器性能曲线绘于图10中,与ES-COD检测器性能进行对比分析。
由图10可以看出:在4种检测器中,ES-COD检测器的检测概率最大;HOC检测器相较于OBF检测器性能有所提高;文献[3]提出的LMS-OBF检测器能够克服海浪磁噪声对航空磁异常检测器的性能影响,相对于HOC检测器,其性能大大提高;在低信噪比条件下,ES-COD检测性能较OBF检测器、HOC检测器和LMS-OBF检测器更优。当信噪比≥-5 dB时,检测概率大于90%,能够准确发现目标,可以满足实际需求;当信噪比≥-7 dB时,检测概率仍保持在50%,相对于文献[3]中的LMS-OBF检测器,ES-COD检测器的弱信号检测能力更突出。基于实测数据分析可知,地磁异常信号与潜艇目标信号相比信号宽度通常较大,飞行器平台干扰信号宽度较小,二者与潜艇目标磁场信号通常呈弱相关性。综上所述,使用联合估计检测方法,在低信噪比条件下能够降低地磁异常和飞行器平台干扰信号引起的虚警概率。
本文针对航空磁探测中潜艇目标信号被强磁背景噪声和干扰淹没导致检测概率低、虚警概率过高的问题,提出了一种航空磁探中的联合估计检测方法,并给出了具体的实现方案。该方法的主要思想是通过磁偶极子建立潜艇目标磁场信号模型,根据观测信号,基于目标先验信息,使用GA和Fréchet距离对目标参量进行估计。基于目标估计信号构造了平均相关检验统计量,设计了使用潜艇目标估计信号的ES-COD联合估计检测器。得到如下结论:
1)使用GA搜索目标参量,能够快速收敛到真实值,得到目标的估计参量。
2)由估计参量得到的潜艇目标磁场估计信号与实测信号相近。
3)实验中ES-COD检测器输出的幅度信噪比约为10,该检测器对航空磁探中的背景磁噪声具有强压制能力,能够显著提高信噪比。
4)当信噪比≥-5 dB时,ES-COD检测器检测概率大于90%,表明该检测器在低信噪比条件下仍具有良好的检测能力。
5)基于奈曼- 皮尔逊准则,在相同虚警概率条件下,ES-COD检测器的检测性能优于OBF检测器、HOC检测器和LMS-OBF检测器,并且可以同时实现目标的检测、定位以及分类,具有一定的实际应用价值。
参考文献(References)
[1] Liu Y, Zhang Y, Yi H. The new magnetic survey method for underwater pipeline detection [J]. Applied Mechanics and Materials, 2013, 239(2): 338-343.
[2] Dames P M, Schwager M, Schwager D. Active magnetic anomaly detection using multiple micro aerial vehicles [J]. IEEE Robotics and Automation Letters, 2016, 1(1): 153-160.
[3] 熊雄, 杨日杰, 王鸿吉. 海浪磁噪声背景中动目标航空磁异常检测算法[J].华中科技大学学报:自然科学版, 2015, 43(5): 100-105.
XIONG Xiong, YANG Ri-jie, WANG Hong-ji. Airborne magnetic anomaly detection algorithm for moving target under ocean wave generated magnetic noise [J]. Journal of Huazhong University of Science and Technology: Natural Science Edition, 2015, 43(5): 100-105. (in Chinese)
[4] Song L, Billings S, Pasion L, et al. Transient electromagnetic scattering of a metallic object buried in underwater sediments [J]. IEEE Transactions on Geoscience and Remote Sensing, 2016, 54(2): 1091-1102.
[5] Liu Z, Pang H, Pan M. Calibration and compensation of geomagnetic vector measurement system and improvement of magnetic anomaly detection [J]. IEEE Transactions on Geoscience and Remote Sensing Letters, 2016, 13(3): 447-451.
[6] Robert F, Saiful H, Mojtaba A, et al. Magnetic signature attenuation of an unmanned aircraft system for aeromagnetic survey [J]. IEEE/ASME Transactions on Mechatronics, 2014, 19(4): 1436-1446.
[7] Ginzburg B, Frumkis L, Kaplan B Z. Processing of magnetic scalar gradiometer signals using orthonormalized functions [J]. Sensors and Actuators, 2002, 102(1): 67-75.
[8] Ginzburg B, Frumkisb L, Kaplan B Z. An efficient method for processing scalar magnetic gradiometer signals [J]. Sensors and Actuators, 2004, 114(1): 73-79.
[9] Frumkisa L, Ginzburg B, Salomonski N. Optimization of scalar magnetic gradiometer signal processing [J]. Sensors and Actuators, 2005, 121(1): 88-94.
[10] Sheinker A, Moldwin M B. Magnetic anomaly detection (MAD) of ferromagnetic pipelines using principal component analysis (PCA) [J]. Measurement Science and Technology, 2016, 27(4): 1-7.
[11] Sheinker A, Salomonskil N, Ginzburg B. Magnetic anomaly detection using entropy filter [J]. Measurement Science and Technology, 2008, 19(5): 1-5.
[12] Sheinker A, Ginzburg B, Salomonski N. Magnetic anomaly detection using high-order crossing method [J]. IEEE Transactions on Geoscience and Remote Sensing, 2012, 50(4): 1095-1103.
[13] Sheinker A, Shkalim A, Salomonski N. Processing of a scalar magnetometer signal contaminated by 1/fαnoise [J]. Sensors and Actuators, 2007, 138(1): 105-111.
[14] Sheinker A, Frumkis L, Ginzburg B. Magnetic anomaly detection using a three-axis magnetometer [J]. IEEE Transactions on Magnetics, 2009, 45(1): 160-167.
[15] 张坚, 林春生, 邓鹏, 等. 非高斯背景噪声下的微弱磁异常信号检测算法[J]. 海军工程大学学报, 2011, 23(4): 22-26.
ZHANG Jian, LIN Chun-sheng,DENG Peng, et al. Detection of weak magnetic anomaly signal in non-Gaussian noise background [J]. Journal of Naval University of Engineering, 2011, 23(4): 22-26. (in Chinese)
[16] 邓鹏, 张坚, 林春生. 航空磁探中水下目标的自适应探测方法[J]. 舰船科学技术, 2012, 34(3): 76-79.
DENG Peng, ZHANG Jian, LIN Chun-sheng. Adaptive detection method of magnetic anomaly detection in aeromagnetic survey [J]. Ship Science and Technology, 2012, 34(3): 76-79. (in Chinese)
[17] 张朝阳, 肖昌汉, 高俊吉, 等. 磁性物体磁偶极子模型适用性的试验研究[J]. 应用基础与工程科学学报, 2010, 18(5): 862-868.
ZHANG Zhao-yang, XIAO Chang-han, GAO Jun-ji, et al. Experiment research of magnetic dipole model applicability for a magnetic object [J]. Journal of Basic Science and Engineering, 2010, 18(5): 862-868. (in Chinese)
[18] 郭志馗, 陈超, 陶春辉, 等. 有限长圆柱体磁异常场全空间正演方法[J]. 地球物理学报, 2017, 60(4): 1557-1570.
GUO Zhi-kui, CHEN Chao, TAO Chun-hui, et al. Forward modeling of magnetic anomalies of finite-length cylinder in 3D space [J]. Chinese Journal of Geophysics, 2017, 60(4): 1557-1570. (in Chinese)
[19] 周家新, 陈建勇, 单志超, 等. 航空磁探中潜艇磁场建模方法分析[J]. 海军航空工程学院学报, 2017, 32(1): 143-148.
ZHOU Jia-xin, CHEN Jian-yong, SHAN Zhi-chao, et al. Analysis of submarine magnetic field modeling method for aeromagnetic detection [J]. Journal of Naval Aeronautical and Astronautical University, 2017, 32(1): 143-148. (in Chinese)
[20] 杨日杰, 熊雄, 郭新奇, 等. 基于潜艇磁偶极子模型的航空磁探潜探测宽度模型与仿真[J]. 兵工学报, 2014, 35(9): 1458-1465.
YANG Ri-jie, XIONG Xiong, GUO Xin-qi, et al. Research on model and simulation of airborne magnetic anomaly detection sweep width based on magnetic dipole model [J]. Acta Armamentarii, 2014, 35(9): 1458-1465. (in Chinese)
[21] Sheinker A, Salomonski N, Ginzburg B, et al. Aeromagnetic search using genetic algorithm [C]∥Progress in Electromagnetics Research Symposium (PIERS 2005). Hangzhou: The Electromagnetics Academy, 2005: 492-495.
[22] 刘渊, 杨永辉, 张春瑞, 等. 一种基于遗传算法的Fuzzing测试用例生成新方法[J]. 电子学报, 2017, 45(3): 552-556.
LIU Yuan, YANG Yong-hui, ZHANG Chun-rui, et al. A novel method for fuzzing test cases generating based on genetic algorithm [J]. Acta Electronica Sinica, 2017, 45(3): 552-556. (in Chinese)
[23] 孙跃, 夏金凤, 唐春森, 等. 采用元胞遗传算法的无线电能传输网路径寻优[J]. 西安交通大学学报, 2017, 51(4): 30-36.
SUN Yue, XIA Jin-feng, TANG Chun-sen, et al. Path optimization for wireless power transfer networks with cellular genetic algorithm [J]. Journal of Xi’an Jiaotong University, 2017, 51(4): 30-36. (in Chinese)
[24] 高翔, 严胜刚, 李斌. 一种混合算法下单个磁通门定位运动磁性目标研究[J].大连理工大学学报, 2016, 56(3): 292-298.
GAO Xiang, YAN Sheng-gang, LI Bin. Study of a hybrid algorithm for localization of mobile magnetic target by a single fluxgate [J]. Journal of Dalian University of Technology, 2016, 56(3): 292-298. (in Chinese)
[25] 高孝杰, 简季, 戴晓爱, 等. 基于Fréchet距离的光谱曲线匹配应用分析[J]. 武汉大学学报:信息科学版, 2016, 41(3): 408-414.
GAO Xiao-jie, JIAN Ji, DAI Xiao-ai, et al. Spectral curve matching application analysis based on Fréchet distance [J]. Geoma-tics and Information Science of Wuhan University, 2016, 41(3): 408-414. (in Chinese)
[26] 安晓亚, 刘平芝, 杨云, 等. 一种线状要素儿何相似性度量方法及其应用[J]. 武汉大学学报:信息科学版, 2015, 40(9): 1225-1229.
AN Xiao-ya, LIU Ping-zhi, YANG Yun, et al. A geometric similarity measurement method and applications to linear feature [J]. Geomatics and Information Science of Wuhan University, 2015, 40(9): 1225-1229. (in Chinese)