在传统的声呐信号处理中,高斯模型占据着主导地位。随着现代主动声呐体制向低频、大功率、大孔径方向发展,混响对声呐性能的影响日趋严重,如何有效地抑制混响干扰成为主动声呐信号处理的首要任务。大量实践表明,沿袭传统的高斯分布假设,往往不能取得理想的抗混响效果。因此,放弃高斯假设而采用非高斯分布假设在水声信号处理业界已成为一种趋势[1-2]。近年来,α稳定分布作为一种非高斯分布模型,具有统计分布的稳定性和概率密度函数(probability density function,PDF)的代数拖尾特点,在信号处理领域受到了广泛的重视,己成为常用的冲激信号统计模型[3-7]。作为唯一的一类构成独立同分布随机变量之和的极限分布[4],α稳定分布是广义的高斯分布,通过调整特征指数α的值实现对PDF拖尾厚度的控制,可以对PDF为单峰钟形随机信号统计特性进行灵活地描述,比高斯分布具有更广泛的适用性。而以混响为主的主动声呐背景干扰,其瞬时幅度样本是均值为零、正负样本量大致相等、上下包络线基本对称的随机信号[1-2],即其PDF恰好是以0为中心左右对称的单峰钟形曲线。因此,本文将原点位置对称α稳定分布(SαS, symmetric alpha-stable distribution)引入到混响非高斯PDF建模中,并使用仿真和实测数据对这一建模方法的可行性和正确性进行检验,最后通过大量海试混响数据建模结果分析,对混响数据的SαS分布规律进行归纳总结。
1 混响的特点与SαS分布简介 1.1 混响数据的特点海洋混响是海洋信道中一切散射作用的总和。从瞬时幅度上看,它是紧跟在发射之后、总体由强至弱的一连串衰减(可能还有起伏震荡)波形。
图 1给出了一段典型的主动声呐阵元级接收数据波形,从中可以分析归纳混响数据的一般特点。
![]() |
Download:
|
图 1 一段典型的主动声呐阵元级接收数据波形 Fig. 1 A typical waveform of active sonar′s receiving data |
引发这段数据的发射信号是4 s脉宽的HFM脉冲。前4 s为发射脉冲直达声,称为直达声区。事实上,发射伊始即引发混响,亦即,直达声区中也有不完整发射脉冲引发的混响,但因其混叠于直达声中难以分离,且此处视为主动声呐近程盲区,一般不对其内的混响进行研究。直达声结束后,数据幅度从与直达声等高度逐渐衰减,直至不再明显降低,大致如图中的4~12 s段所示,这一段是典型混响区域,称为混响限制区,该区域内,混响强度高于背景噪声,是声呐的主要背景干扰。大致从12 s向后,混响淹没于背景噪声之中,以海洋环境噪声为主的背景噪声成为主要的背景干扰,混响的影响可以忽略,而背景噪声是宽平稳的,故接收数据幅度不再明显变化,此段称为噪声限制区。根据声呐量程,噪声限制区可以延续很长,这里仅给出一部分。不难理解,研究混响数据的统计特性,要基于混响限制区内的数据进行。
从图 1中还直观可见,混响数据是0均值的、正负样本量大致相等、上下包络线基本对称,这意味着,可以使用0均值对称单钟型曲线来对混响数据的概率密度函数进行拟合建模。零均值高斯分布、混合高斯分布,和下面要介绍的SαS分布,其PDF恰好就是这样一类曲线。
1.2 SαS分布简介若随机变量X具有如下形式的特征函数[3]:
$ \begin{array}{*{20}{c}} {\varphi (t) = E[\exp ({\rm{i}}tX)] = }\\ {\exp \left\{ {{\rm{i}}\mu t - \gamma |t{|^\alpha }[1 + {\rm{i}}\beta {\mathop{\rm sgn}} (t)\omega (t,\alpha )]} \right\}} \end{array} $ | (1) |
则称X服从α稳定分布,记为
$ \omega (t,\alpha ) = \left\{ {\begin{array}{*{20}{l}} {\tan \left( {\frac{{{\rm{ \mathsf{ π} }}\alpha }}{2}} \right),\alpha \ne 1}\\ { - \frac{2}{{\rm{ \mathsf{ π} }}}\ln |t|,\alpha = 1} \end{array},{\mathop{\rm sgn}} (t) = \left\{ {\begin{array}{*{20}{c}} {1,}&{t > 0}\\ {0,}&{t = 0}\\ { - 1,}&{t < 0} \end{array}} \right.} \right. $ |
式中:α为特征指数,表征稳定分布脉冲特性强弱,α∈(0, 2],值越小,则所对应分布的拖尾越厚,非高斯脉冲特性越显著;γ为离差,亦称为尺度参数、分散系数,表征稳定分布随机变量偏离其均值或中值的程度,γ∈(0, ∞),定义σ=γ1/α为标准离差;β为偏斜指数,表这种分布称为原点征稳定分布的偏斜程度,β∈[-1, 1],β=0、β>0和β<0时,分别对应稳定分布无偏斜(对称)、右偏斜和左偏斜情况;μ为位置参数,表征稳定分布的均值(α>1时)或中值(α≤1时),μ∈(-∞, ∞)。
α稳定分布特征函数与其PDF是一对傅立叶变换[3],可通过对特征函数进行傅立叶变换得到PDF:
$ \begin{array}{*{20}{c}} {f(x;\alpha ,\beta ,\gamma ,\mu ) = }\\ {\frac{1}{{\rm{ \mathsf{ π} }}}\int_0^\infty {\exp } \left\{ { - \gamma {t^\alpha }} \right\}\cos \left[ {(x - \mu )t + \gamma {t^\alpha }\beta \omega (t,\alpha )} \right]{\rm{d}}t} \end{array} $ | (2) |
除几组特殊参数取值外,式(2)一般没有封闭的解析表达式,需要以数值积分方式计算[8-9]。图 2给出了以这种方式计算绘出的不同参数取值时的α稳定分布的PDF曲线比较。从各图中参数变化时的PDF曲线对比,可以清晰地看到对应参数的物理意义体现。
![]() |
Download:
|
图 2 不同参数下的S(α, β, γ, μ)分布PDF变化 Fig. 2 PDFs of S(α, β, γ, μ) with different parameters |
一些常见分布可以看作是α稳定分布S(α, β, γ, μ)的特例[3],比如:当α=2、β=0时,S(2, 0, γ, μ)就是高斯分布(μ, 2γ);当α=1、β=0时,S(1, 0, γ, μ)就是柯西分布(μ, γ)。
当偏斜指数β=0时,α稳定分布称为对称α稳定分布,记为SαS分布。不难发现,高斯分布、柯西分布都是SαS分布的特例。正常状态下的主动声呐背景干扰数据,包括混响数据,其样本都是无偏斜且0均值的,故可用位置参数μ=0的SαS分布来拟合其PDF,可把这种分布称为原点位置SαS分布,简记作Sα(γ)。
可以看到,Sα(γ)分布仅有2个参数:特征指数α和离差γ,就像常用的高斯分布一样简洁。值得指出,在对混响数据进行非高斯PDF建模时,堪用的模型,即使最简洁的二阶0均值混合高斯模型ZMGM(2),也有3个参数[1]。相对少的参数意味着相对快捷而精确的参数估计算法,在这一点上,Sα(γ)具有无可比拟的优势。
2 混响数据SαS分布建模方法获取合格而足量的混响样本数据后,即可以原点位置SαS分布对其进行建模。而建模的核心任务是,如何构建Sα(γ)的2个未知参数α和γ估计器。
一般的α稳定分布S(α, β, γ, μ)参数估计问题,已经有很多文献研究[3, 8, 10]。如前分析,混响数据实际上可以简化的α稳定分布,即Sα(γ) 分布进行PDF建模,故针对这一特点,本文从诸多方法中予以选取并简化改进,建立了如下2种方法:正负矩法和对数矩法。
2.1 正负矩法设随机变量X~Sα(γ),则有[3]:
$ \begin{array}{l} {\rm{E}}\left( {|X{|^p}} \right){\rm{E}}\left( {|X{|^{ - p}}} \right) = C(p,\alpha )C( - p,\alpha ) = \\ {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} \frac{{2\tan (p{\rm{ \mathsf{ π} }}/2)}}{{\alpha \sin (p{\rm{ \mathsf{ π} }}/\alpha )}},0 < p < \min (\alpha ,1) \end{array} $ | (3) |
式中:
$ {{M_{|X|}}(p) \buildrel \Delta \over = {\rm{E}}\left( {|X{|^p}} \right) = C(p,\alpha ){\gamma ^{p/\alpha }}} $ | (4) |
$ {C(p,\alpha ) = \frac{{{2^{p + 1}}\varGamma [(p + 1)/2]\varGamma ( - p/\alpha )}}{{\alpha \sqrt {\rm{ \mathsf{ π} }} \varGamma ( - p/2)}}} $ | (5) |
这里,Γ(·)为伽马函数,定义为[8]:
$ \varGamma (x) \buildrel \Delta \over = \int_0^\infty {{t^{x - 1}}} {{\rm{e}}^{ - t}}{\rm{d}}t $ | (6) |
由式(3)~(5),有:
$ \begin{array}{*{20}{c}} {{\mathop{\rm sinc}\nolimits} \left( {\frac{{p{\rm{ \mathsf{ π} }}}}{\alpha }} \right) = \frac{{\sin (p{\rm{ \mathsf{ π} }}/\alpha )}}{{(p{\rm{ \mathsf{ π} }}/\alpha )}} = }\\ {\frac{{2\tan (p{\rm{ \mathsf{ π} }}/2)}}{{p{\rm{ \mathsf{ π} }}{M_{|X|}}(p){M_{|X|}}( - p)}},0 < p < \min (\alpha ,1)} \end{array} $ | (7) |
通过样本统计得到
阶数p的设定是正负矩法的关键之一,但多数文献中仅给出一范围指示[3, 10],并无深入探讨。针对混响建模实际,本文采取如下做法:首先要根据经验确定待建模数据的特征指数下限值αmin(这一经验值的获取可通过下文所述的对数矩法建模总结),确保p取值小于αmin或1;为提高估计精度,可设置多组阶数值pi通过式(7)得到多组估值
设随机变量X~Sα(γ),则由式(4)有:
$ {\rm{E}}\left( {|X{|^p}} \right) = {\rm{E}}\left( {{{\rm{e}}^{p\ln |X|}}} \right) = C(p,\alpha ){\gamma ^{p/\alpha }} $ | (8) |
引入负阶矩和分数阶矩概念后,式(8)在p=0处连续[3]。定义随机变量Y=ln|X|,则Y的矩母函数为:
$ {\rm{E}}\left( {{{\rm{e}}^{pY}}} \right) \buildrel \Delta \over = \sum\limits_{k = 0}^\infty {\rm{E}} \left( {{Y^k}} \right)\frac{{{p^k}}}{{k!}} = C(p,\alpha ){\gamma ^{p/\alpha }} $ | (9) |
故Y的任意k阶矩都是有限的,且满足:
$ {\rm{E}}\left( {{Y^k}} \right) = {\left. {\frac{{{{\rm{d}}^k}}}{{{\rm{d}}{p^k}}}\left[ {C(p,\alpha ){\gamma ^{p/\alpha }}} \right]} \right|_{p = 0}} $ | (10) |
k=1时可得Y 的一阶原点矩(均值)为:
$ {\mu _Y} = {\rm{E}}(Y) = \left( {\frac{1}{\alpha } - 1} \right){c_e} + \frac{1}{\alpha }\ln \gamma $ | (11) |
式中ce为Euler常数。同理可得Y的二阶中心距(方差)为:
$ \sigma _Y^2 = {\rm{E}}\left\{ {{{[Y - {\rm{E}}(Y)]}^2}} \right\} = \frac{{{{\rm{ \mathsf{ π} }}^2}}}{6}\left( {\frac{1}{{{\alpha ^2}}} + \frac{1}{2}} \right) $ | (12) |
混响建模实际应用中,可首先根据样本统计得到Y=ln|X|的均值和方差估计
应当指出,式(11)与(12)中已无矩母函数自变量(亦即正负矩法的中阶数)p,这就是说,对数矩法不需要p初值设定,恰与正负矩法形成补充,二者估计结果可以相互印证。
2.3 仿真实例首先使用PDF可设定的仿真数据检验Sα(γ) 分布对高斯和非高斯PDF的拟合能力。
图 3是一段标准正态分布随机序列的建模拟合情况。其中图 3(a)为波形图,图 3(b)为PDF曲线对比,虚线为使用设定的高斯分布参数[μ, σ2] =[0, 1]绘出的曲线,代表PDF真值;实线为使用对数矩法Sα(γ)分布参数估计[
![]() |
Download:
|
图 3 SαS分布Sα(γ)对高斯数据PDF的拟合性能 Fig. 3 PDF fitting performance of Sα(γ) for Gaussian data |
图 4是一段非高斯随机序列的建模拟合情况。非高斯序列按参数[εI, σB2, σI2] =[0.1, 1, 100]的混合高斯分布以复合抽样法[1]生成。(a)为波形图,从中可观察到大量“毛刺”(尖峰脉冲)突出,故这是一种比较严重的非高斯情形。图 4(b)为PDF曲线对比,虚线为使用设定的混合高斯参数绘出的曲线,代表PDF真值;细实线为使用对数矩法Sα(γ)分布参数估计(
![]() |
Download:
|
图 4 SαS分布Sα(γ)对非高斯数据PDF的拟合性能 Fig. 4 PDF fitting performance of Sα(γ) for non-Gaussian data |
本节以数组低频主动拖曳线列阵声呐(0411三亚海试)数据分析为例,验证以混响为主的主动声呐背景干扰数据(即图 1所示中的混响限制区内数据,简称为混响)SαS分布PDF建模方法的正确性,并初步归纳混响PDF的SαS分布规律。
3.1 建模方法验证(0411三亚海试)发射3种主动脉冲信号:Ⅰ表示脉宽为432 ms、频率为750 Hz的CW脉冲,称为窄带短脉冲;Ⅱ表示同样脉宽、频带为650~850 Hz的LFM脉冲,称为宽带短脉冲;Ⅲ表示脉宽为1 296 ms、同样频带的LFM脉冲,称为宽带长脉冲。使用80元拖曳线列阵接收回波数据。进行PDF建模前数据已经过带通滤波→波束形成→移频→降采样→幅值归一化等处理[1, 11]。
图 5~7分别给出了3段3种脉冲混响数据的波形图和SαS分布PDF建模拟合图。为便于比较,图 5~7(b)图中同时绘出了:①使用柱状图统计法得到的PDF曲线,图中标记为HS。在足够样本量下,HS曲线可视为真实分布。②使用正负矩法Sα(γ) 参数估计得到PDF拟合曲线,图中标记为NM。③使用对数矩法Sα(γ) 参数估计得到PDF拟合曲线,图中标记为LM。NM和LM方法得到的α、γ参数估计值也同时标记于图 5~7(b)图上方。
![]() |
Download:
|
图 5 一段窄带短脉冲混响数据SαS分布验模 Fig. 5 SαS PDF fitting performance for a reverberation segment of narrow-band short pulse |
![]() |
Download:
|
图 6 一段宽带短脉冲混响数据SαS分布验模 Fig. 6 SαS PDF fitting performance for a reverberation segment of wide-band short pulse |
![]() |
Download:
|
图 7 一段宽带长脉冲数据SαS分布验模 Fig. 7 SαS PDF fitting performance for a reverberation segment of wide-band long pulse |
由图可见,负阶矩法、对数矩法2种参数估计方法得到的参数估计值几乎完全相同,拟合的NM和LM曲线几乎完全重合。而此二者又与统计得到的HS曲线几乎完全重合。这说明,这2种参数估计方法具有几乎完全相同的PDF建模精度,都能够真实反映数据的PDF特性。换言之,SαS分布Sα(γ) 可以精确反映混响数据的PDF,其参数估计选用正负矩法或对数矩法皆可。从运算速度看,对数矩法略优。下文选用对数矩法。
3.2 建模结果分析选取(0411三亚海试)中不同方位通道、不同距离处的Ⅰ、Ⅱ、Ⅲ 3种发射信号引起的混响数据各5 430段,以对数矩法进行SαS分布建模。由于离差估值
表 1给出了这些数据的建模结果。3组共计16 290段数据的
![]() |
表 1 海试混响数据α估值分布区间统计 Table 1 Statistics of α-estimation intervals for the sea trial data |
由表可见,各组数据的
因此,可以说,混响数据的SαS分布特征指数近半数位于1.6~1.9,表现出一定的非高斯性。相较而言,宽带长脉冲信号引发的混响数据非高斯性略强,而窄带短脉冲信号引发的混响数据的特征指数分布范围更广(信号Ⅰ组有少量
1)SαS分布可对主动声呐背景干扰进行较好的概率密度拟合;
2) 在混响限制区内,主动声呐背景干扰的特征指数α近半数集中在1.6~1.9,表现出一定的非高斯性。
[1] |
王平波. 主动声纳非高斯信号处理技术研究[D]. 武汉: 海军工程大学, 2006. WANG Pingbo. Study on non-Gaussian signal processing in active sonar[D]. Wuhan: Naval University of Engineering, 2006. ( ![]() |
[2] |
王汗青. 基于对称Alpha稳定分布的水声信号处理技术研究[D]. 武汉: 海军工程大学, 2012. WANG Hanqing. Study on underwater acoustic signals processing based on symmetrical Alpha stable distribution[D]. Wuhan: Naval University of Engineering, 2012. ( ![]() |
[3] |
邱天爽, 张旭秀, 李小兵, 等. 统计信号处理——非高斯信号处理及其应用[M]. 北京: 电子工业出版社, 2004. QIU Tianshuang, ZHANG Xuxiu, LI Xiaobing, et al. Statistical signal processing: non-Gaussian signal processing and application[M]. Beijing: Publishing House of Electronics Industry, 2004. ( ![]() |
[4] |
SHAO M, NIKIAS C L. Signal processing with fractional lower order moments: stable processes and their applications[J]. Proceedings of the IEEE, 1993, 81(7): 986-1010. DOI:10.1109/5.231338 ( ![]() |
[5] |
SHAO Min, NIKIAS C L. Signal detection in impulsive noise based on stable distributions[C]//Proceedings of 27th Asilomar Conference on Signals, Systems and Computers. Pacific Grove, CA, USA, 1993: 218-222.
( ![]() |
[6] |
RUPI M, TSAKALIDES P, NIKIAS C L, et al. Robust constant modulus arrays based on fractional lower-order statistics[C]//1999 IEEE International Conference on Acoustics, Speech, and Signal Processing. Phoenix, AZ, USA, 1999: 2945-2948.
( ![]() |
[7] |
陈思佳. Alpha稳定分布噪声下的自适应滤波算法研究[D]. 杭州: 杭州电子科技大学, 2020. CHEN Sijia. Research on adaptive filtering algorithm under alpha stable distribution noise[D]. Hangzhou: Hangzhou Dianzi University, 2020. ( ![]() |
[8] |
彭成. 基于分数低阶统计量的非高信号处理技术研究[D]. 武汉: 海军工程大学, 2014. PENG Cheng. Study on non-Gaussian signal processing based on fractional lower order statistics[D]. Wuhan: Naval University of Engineering, 2014. ( ![]() |
[9] |
彭成, 王平波, 刘旺锁. Alpha稳定分布序列的数值仿真方法[J]. 声学技术, 2014, 33(5): 473-476. PENG Cheng, WANG Pingbo, LIU Wangsuo. Numerical simulation method of Alpha stable distribution sequence[J]. Technical acoustics, 2014, 33(5): 473-476. ( ![]() |
[10] |
王余, 王本猛, 卫红凯, 等. 对称α稳定分布特性及其参数估计方法研究[J]. 舰船电子工程, 2017, 37(4): 116-121. WANG Yu, WANG Benmeng, WEI Hongkai, et al. Symmetrical Alpha stable distribution and its parameter estimation method[J]. Ship electronic engineering, 2017, 37(4): 116-121. ( ![]() |
[11] |
WANG Pingbo, WEI Hongkai, HU Jinhua, et al. An adaptive beamforming algorithm for underwater acoustic signal based on α-stable distribution[C]//2014 IEEE Workshop on Electronics, Computer and Applications. Ottawa, ON, Canada, 2014: 657-659.
( ![]() |