脉冲激光雷达信号降噪方法对比

陈冬，王江安，康圣

（海军工程大学电子工程学院光电研究所，湖北武汉 430033）

0 引言

在激光遥感领域，研究者的注意力主要集中在信号反演方法上而忽视了降噪处理算法的深入研究。通常采用的仅仅是中值滤波、窄带滤波等传统方法，也有研究者采用多次测量取均值法。近年才有研究者应用小波阈值法进行滤波并取得了较好效果［1－2］。

对脉冲激光雷达信号频谱进行分析可知，其近场信号频谱成分分布广泛而远场高频分量能量较低而低频分量起主导作用。因而较之传统的滤波方法，时频分析工具等现代信号处理方法更适合于该课题。经验模态分解理论是美籍华人黄鄂教授提出的一种基于时频分析的新算法，它与小波分解算法的区别在于对瞬时频率定义的差别。本文采用经验模态分解（EMD）滤波对激光雷达信号进行降噪处理，并与扩展卡尔曼算法（EKF）和小波滤波进行比较。将3种方法在不同大气条件下的处理结果进行对比，从而得到最优的降噪方法。

1 信号模型与信号质量评价

1.1 信号模型

仿真验证时如果采用脉冲激光雷达回波实测数据，无法对原始数据和滤波后数据的信号进行评价，因而必须构建仿真理想信号。严格地说，每种大气条件下的信号模型都是不同的。它受到大气粒子尺寸分布，粒子形状分布、浓度，粒子性质及激光雷达发射和接收视场角等多种因素影响。研究者通常把大气分为3种情况进行建模，把非常稀薄和非常稠密的粒子群分布单独考虑，分别采用单粒子散射近似和漫反射近似，而最普遍的情况是采用一级多次散射近似进行建模。大气脉冲激光雷达方程（式（1））正是在一级多次散射近似条件下推导出来的，适用于绝大多数大气回波情况，也与实测信号最为接近。

式中:Pr（R）为距离R处的回波功率;βπ（R）为后向散射系数;σ（R）为消光系数;Cq为仪器常数，它又可以表示为:

式中:Ar为接收镜头面积;Tt和Tr分别为发射/接收光学系统透过率;E0为发射脉冲能量;τp为脉冲宽度;c为光速;η（R）为几何重叠因子。

后向散射系数是激光雷达方程中的重要参数，根据Klett的研究［3］，在米散射起主导作用时，消光系数与后向散射系数的关系为β（z）=B0σk（z），B0为后向散射比，k的取值和波长与粒子尺度的相对值有关，通常取值区间为［0．66，1］。

实际的激光雷达接收信号要通过光电探测器、高频放大器、模数转换器等器件。根据实际器件参数，将信号模型进一步修正。设放大器的偏置电平为lv，AD转换器的采样频率为f，同时假定第n－1个采样点与第n个采样点范围内采样体积内的平均消光系数为，回波模型可以修正为:

图1 仿真、实测回波信号比较Fig.1Compare between real and simulated signal

为了检验模型的正确性，将真实信号与同等条件仿真波形进行比较。图1中虚线为某型激光雷达在某机场晴朗均匀大气中的实测回波信号，实线是采用式（3）拟合的理想回波信号。可以看出，仿真信号对远场回波信号拟合较好，近场信号存在差别，这是几何重叠因子的差异造成的。仪器常数的诸多因子中，几何重叠因子是难以确定的。有研究人员假设光束能量高斯分布通过激光光束发射与接收视场截面的大小来确定几何重叠因子的大小［4］。但激光光束能量的分布特性的随机性大，近场能量强、多次散射强烈使得这种方法意义不大。因此，研究中只关注视场完全重叠后的区域，近场拟合的精确性对仿真效果并无影响。

1.2 噪声模型

激光雷达的噪声按形成原因可分为散粒噪声、热噪声、背景光噪声及1/f噪声等［2］。其中除了1/f噪声外，其他噪声都可认为是加性噪声，并且符合高斯分布。因此研究中采用高斯噪声模型。

1.3 信号质量评价方法

评价信号质量通常采用的参数是信噪比，但信号信噪比的定义有很多，本文选用信号幅值均方根与噪声幅度均方根之比来定义。

式中:Pi为信号幅度;Ni为噪声幅度;n为信号长度。

2 降噪算法

2.1 EKF滤波

在激光雷达领域，研究者将扩展卡尔曼算法当作1种大气状态参数的反演方法［5－6］。这样难以对算法降噪性能加以评估，也无法验证反演模型的准确性。因此笔者倾向于把EKF法修正为1种降噪算法。

将消光系数作为信号主要参量，信号模型为

其中:v（n）为消光系数波动参数，它服从正态分布，方差为R，均值为0。

偏置电平很容易滤除，去掉放大器偏置电平并考虑加性噪声，式（3）可化简为式（6），其中w（n）均值为0，方差为Q。

初始条件为P（n0）=0，n0为发射与接收视场的完全重叠位置，此时η（n0）=1。

2.2 小波滤波

小波滤波采用多尺度正交小波基对信号进行分解，将信号和噪声在不同时频域进行分离，然后在各时频范围内对依据阈值进行滤波，最后重构信号从而达到降噪的目的。离散小波基的定义为:

任何函数都可以分解为小波基的加权和，

其中，加权值即为小波系数:

小波理论产生于20世纪80年代，学者们对小波降噪方法的研究广泛且深入，因此本文不再讨论。对于激光雷达信号降噪问题，小波基和阈值的选取是最重要的2个方面。不同小波基的特性不同，选取不同的小波基对滤波结果影响很大。Db5小波具有双正交性、紧支撑性和近似对称性，十分适合对指数衰减型曲线进行滤波。阈值的选取对滤波效果影响也很大，常用的降噪阈值选取方法有长度对数阈值（Sqtwolog）、最小极大方差阈值（Minimaxi）、无偏风险阈值（Rigrsure）、启发式阈值（Heursure）等［8－10］。这些方法的主要思想是根据小波基长度与噪声均方根构造最佳阈值，它们都可以改善信号信噪比，但降噪能力有限。对于本课题来说，仿真表明固定阈值降噪效果更好，当阈值选择为二倍噪声均方根时可使滤波效果达到最佳。

2．3 EMD分解滤波

EMD分解滤波也是一种基于时频分析的算法，它与小波分解算法的区别在于对瞬时频率定义的差别。EMD算法具有自适应性，它会为每条信号曲线产生不同的正交基。

经验模态分解（EMD），是通过固有模态函数（IMF）提取来实现的［11］。1个函数可以分解为n个固有模态函数ci（t）与1个剩余分量r（t）的和，即

固有模态函数IMF必须满足以下2个条件:一是信号持续时间内，其过0点的个数和极值点的个数相等或最多相差1个;二是信号极大值或极小值构成的包络线的均值为0。

固有模态函数采用以下步骤进行提取:寻找信号f（t）的局部极值点;在极大值之间进行插值拟合处理形成最大值包络曲线maxi（t），在极小值之间进行插值拟合处理形成最小值包络曲线mini（t）;求包络均值曲线m（t）=［max（t）+min（t）］/2;信号减去包络均值曲线;多次迭代后，剩余信号满足停止准则时筛选结束。筛选门限见式（15），筛选门限一般在0.2～0.3之间。

经过EMD提取之后，信号被分解为多个固有模态函数，对固有函数进行阈值滤波，将滤波后的数据重构就可以达到降噪的目的。

仿真表明，筛选门限取值2.5时滤波效果较好，软阈值法的效果比硬阈值法滤波效果好，当阈值为噪声均方根时可使降噪效果达到最佳。

3 仿真对比

仿真主要针对2种大气条件下的回波信号:均匀大气和局部突变大气。均匀大气主要针对无风或微风条件下的自然大气状态，这时消光系数比较稳定随距离起伏不大。局部突变大气针对存在风切变、湍流及云雾或烟尘等天气条件，它们会改变大气的局部分布，造成局部消光系数突变。

3.1 均匀大气条件

常规大气条件下，消光系数均匀，在均值附近上下起伏。因而选择不同的消光系数均值和方差构造激光脉冲回波，并叠加高斯噪声。如图2和图3所示。

仿真依据消光系数从低至高进行，每个消光系数取值又采用不同的噪声幅度进行仿真，仿真结果如表2所示。

表2中第1行为仿真消光系数值，第1列为噪声均方根（δ）。表中数据依次为原始信号EKF滤波小波滤波EMD滤波对应的信噪比。

3.2 局部突变条件

此种条件时，绝大部分消光系数分布与均匀大气下相同，只是局部空气中消光系数陡增骤减。当消光系数如图4所示分布时，回波信号波形如图5所示，仿真结果如表3所示。

表3中第1行为仿真消光系数值，第1列为噪声均方根（δ）。表中数据依次为原始信号EKF滤波小波滤波EMD滤波对应的信噪比。

4 结语

扩展卡尔曼滤波的降噪能力有限，但可以将信号直接反演为消光系数，初始值的选取对这种方法的滤波效果影响很大。消光系数较高时，EKF滤波后的失真较大;在消光系数较低时的滤波效果较好。在均匀大气条件下，EKF法反演精度较高;而当大气中有局部突变时，若消光系数陡增幅度过大会破坏这种方法的成立条件，导致误差较大。小波滤波的效果在多数情况下是几种方法中的最佳方法，信噪比较高时，小波滤波的优势并不明显，并且有可能造成信号失真;中低信噪比条件下，这种方法对信号质量的改善最为明显;在突变大气条件下，尽管这种方法的滤波效果降低了，但仍然是最佳方法。在信噪比较高时EMD滤波效果明显好于其他2种方法。目前采用EMD法对激光雷达回波信号进行滤波的研究仍处于起步阶段，在降噪阈值和筛选门限这2个方面进一步研究会使该方法的优势得到进一步体现。

综上所述，本文根据激光雷达性能参数和典型的大气条件建立了1种更为实用的回波信号仿真模型。本文工作证明EMD滤波法在某些大气条件下可以达到最佳滤波效果，是现有方法的必要补充。综合采用3种方法才能使滤波达到最佳效果。均匀大气条件下，如果噪声幅度较高，应采用小波滤波;噪声幅度较低时，如果消光系数较小，应采用EKF滤波;如果消光系数较大，就应采用EMD算法。突变大气条件下，如果噪声幅度较高，宜采用小波滤波;如果噪声幅度较低则应采用EMD滤波。将该方法应用到激光雷达消光系数测量工程中，得到了较好的效果。

［1］FANG Hai-tao，HUANG De-shuang．Noise reduction in lidar signal based on discrete wavelet transform［J］．Optic communication，2004，233:67－76．

［2］Aime Lay-Ekuakille，TROTTAA．Comparisonbetween adjustablewindowtechniqueandwaveletmethodin processing backscattering lidar sinal［C］．SPIE，2004 （5240）:72－82．

［3］KLETT J D．Stableanalyticalinversionsolutionfor processing lidar returns［J］．Applied Optics．1981，20（2）: 211－210．

［4］王治华，贺应红，左浩毅，等．基于高斯光束特性的Mie散射大气激光雷达回波近场信号校正研究［J］．物理学报，2006，55（6）:3188－3192．

WAN Zhi-hua，HE Ying-hong，ZUO Hao-yi，et al．The correction of short-range Mie scattering laser lidar returns based on the Gaussian character of laser beam［J］．Acta Physica Sinica，2006，55（6）:3188－3192．

［5］LAINIOTIS D G，PAPAPARASKEVA P，KOTHAPALLI G，et al．Adaptive filter applications to LIDAR:return power and log power estimation［J］．IEEE Transactons on Geoscience and Remote Sensing，1996，34（4）:886－891．

［6］José M．Bioucas Dias，et al．Reconstruction of backscatter and extinction coeffictents in lidar:A stochastic filtering approach［J］．IEEE Transactions on Geoscience and Remote Sensing，2004，42（2）:443－456．

［7］HUANG N E，ZHENG Shen，STEVEN R L，et al．The empirical mode decomposition and the Hilbert spectrum for nonlinear and non-stationary time series analysis［J］．The Royal Society London A，1998，454:903－995．

［8］张斌，王彤，谷传纲，等．改进的小波阈值消噪法在湍流信号处理中的应用［J］．工程热物理学报，2009，30（3）: 407－410．

ZHANG Bin，WANG Tong，GU Chuan-gang，et al．Applcation of improved wavelet threshold de-noising in turbulent singnal processing［J］．Journal of Engineering Thermophysics，2009，30（3）:407－410．

［9］DONOHO D L．Denoising by soft-thresholding［J］．IEEE Trans on Information Theory，1995，41（3）:613－627．

［10］DONOHO D L，JOHNSTONE I M．Ideal spatial adaptation by wavelet shrinkage［J］．Journal of the Royal Statistical Society，1997，59（2）:319－351．