时频分析在激光雷达中的应用进展 下载: 886次
1 引言
激光雷达具有精度高、时空分辨率高、抗干扰能力强、探测范围广及可多目标探测等突出优点,在遥感探测、大气参数探测、航空航天、海洋和**等领域具有重要的应用[1-5]。经典的傅里叶变换从独立域出发对信号进行描述,不能刻画信号频率随时间变化的情况,适用于分析固定信号。而激光雷达回波信号是典型的非线性非平稳信号,时频分析正是为了描述信号时变频率分量而发展起来的一种分析方法。相比于传统的时域或频域处理,时频分析从时间-频率域揭示信号全貌,通过构造时间和频率的联合密度函数,将一维域分析推广至二维(2D)时频空间进行描述,对处理非线性非平稳信号具有重要意义[6]。时频分析在激光雷达的大气参数特征分析与提取、信号去噪、运动目标成像与检测、微多普勒特征分析等方面已得到广泛研究。
时频分析主要分为线性时频分布、双线性时频分布、自适应非参数类时频分布和自适应参数类时频分布等[7-8]。
线性时频分布是由傅里叶变换演化而来,满足线性叠加性,不存在交叉项。常见的线性时频表示主要有Gabor展开[9]、短时傅里叶变换[10]、小波变换[11]、S变换[12]、分数阶傅里叶变换[13]等。Gabor展开是最早提出的信号时频联合表示的一种方式,Gabor展开的基本原理是用窗函数及其时移和频移的函数组进行展开,展开系数就是Gabor变换,Gabor展开对信号采用时间和频率的二维描述方式,被称为“信息图”。短时傅里叶变换的基本思想是利用一个固定大小的可滑动时间窗将原始信号截断为许多小的时间间隔,截断的窗内信号近似平稳,作每段信号的傅里叶变换,得到各时段的信号频谱,故又称为加窗傅里叶变换。Gabor变换和短时傅里叶变换采用固定窗函数,存在时间分辨率和频率分辨率之间的权衡。与短时傅里叶变换采用单一的分析窗不同,小波变换在低频处使用长窗函数,在高频处使用短窗函数,具有可变分辨率,在多尺度细化分析方面优于短时傅里叶变换。S变换作为短时傅里叶变换和小波分析的继承和发展,克服了短时傅里叶变换单一窗函数的缺点,又具有小波变换多尺度分析的特性。S变换采用窗宽与频率的倒数成正比的高斯窗函数,低频时具有较高的频率分辨率和较低的时间分辨率,高频时具有较高的时间分辨率和较低的频率分辨率。分数傅里叶变换是一种广义的时频分析,在时频平面内,分数傅里叶变换是一种坐标轴绕原点逆时针旋转任意角度后构成的分数阶傅里叶域的表示方法。
双线性时频分布也被称为二次型时频分布,反映的是信号在时频面的能量分布,它是信号的二次变换,不满足线性叠加性。双线性时频分布包括Cohen类时频分布[14]、Affine类时频分布[15]、重排类时频分布[16]和自适应最优核时频分布[17]。Wigner-Ville分布是分析非平稳信号最具代表性的双线性时频分布,对线性调频信号是无偏估计,具有最优分辨率,同时满足许多理想的数学性质如能量守恒、实值、对称性、时-频边缘特性、时移和频移不变性、时频伸缩和有限支撑特性等。但Wigner-Ville分布会使多分量信号和非线性信号产生严重交叉项,严重干扰信号时变谱的可分辨性和可解释性[18-19]。围绕着减少或消除Wigner-Ville分布中的交叉项,Cohen[14]对已有的时频方法研究成果进行总结,发现大多数的双线性时频分布都可以一个普遍的形式来表示,不同的核函数对应不同的时频分布。常用的Cohen类时频分布有谱图[20]、Wigner-Ville分布、伪Wigner-Ville分布[21]、平滑Wigner-Ville分布[22]、Born-Jordan分布[23]、Choi-Williams分布[24]、Zhao-Atlas-Marks分布[25]、Page分布[26]、Rihaczek分布[27]、Margenau-Hill分布[28]等。Affine类时频分布[15]中最著名的分布当推尺度图,Wigner-Ville分布为连接Cohen类和Affine类时频分布的桥梁,前者基于Wigner-Ville分布进行时频平滑,后者基于时间-尺度平滑处理。重排类时频分布[16]的思想认为本地化能量类比于质量分布,无需对称分布在几何中心处,应分配在更能代表本地化能量分布的重心处。重排方法有效提高了信号分量的时频聚集性,同时减小交叉项,改善了时频面的可解释性。自适应最优核时频分布[17]采用短时模糊函数和随时间而自适应改变的核函数,使自适应项最大程度地被保留,同时使多分量信号中的交叉项最大程度地被抑制。
自适应非参数类时频分布是随信号特征自适应改变的时频分布,主要代表有经验模式分解(EMD)[29]以及全新的信号处理方法Hilbert-Huang变换[30]。经验模式分解可以自适应地从信号中得到基函数,根据信号尺度不同,将信号分解为从瞬态尺度到粗糙尺度的不同的信号分量,实际上信号中的噪声往往分布在高频区域,常采用舍弃高频部分来去除噪声。Hilbert-Huang变换通过经验模式分解将任意复杂的信号分解为若干个本征模态(IMF)之和,对每个IMF进行Hilbert变换获得具有物理意义的瞬时频率和振幅,从而得到非平稳信号的时频分布。
自适应参数类时频分布采用时频基函数,自适应地选择基函数参数将信号分解到联合时频域。自适应参数化时频分布主要有自适应高斯表示法[31]、匹配追踪算法[32]、自适应chirplet分解[33]等。自适应高斯表示法将时间分辨率、频率分辨率和时频中心调整到与信号最匹配的程度,采用归一化高斯函数来展开信号,归一化高斯函数的标准差是可调节的。匹配追踪算法用不同波形特征的时频原子构造庞大的“字典”,用匹配追踪算法从中挑选与信号波形匹配的原子,从而达到从信号中提取信息的目的。自适应chirplet分解将具有时间中心、频率中心和时频范围的chirp函数作为基函数,搜索到信号基函数组中所有可能基上的最大投影并自适应选择它们的参数。
2 时频分析在激光雷达中的应用
2.1 大气参数特征分析与提取
激光雷达的基本原理是出射激光脉冲与大气相互作用,采用光学天线收集大气后向散射信号,经光学接收机处理、探测和数据收集后,根据时频分析算法可以反演得到一系列关键大气参数,诸如风速、重力波波动、云高、PM2.5、气溶胶浓度、温度、湿度、能见度及大气成分等。通过时频分析可以提取出信号的瞬时频率曲线、矩、边缘特性、Renyi信息、小波系数、小波包能量、Hilbert谱及边际谱等信息。
2010年,鲁汶大学采用双线性时频算法来探测气象灾害,如飞机尾部涡流和风切变等[34]。考虑了大气湍流效应和地面效应,对激光在湍流介质中的传播进行数值模拟仿真。该研究由Cohen类时频算子导出相干激光雷达Wigner-Ville谱的分析方程,并推导用于探测飞机尾涡的自适应谱模型,该模型可以用来估计风速似然密度方程和极值径向速度廓线。
图 1. 尾涡仿真结果和径向廓线的Wigner-Ville分布。(a)尾涡对中径向风速的等高线数值模拟图;(b)三条视线的径向风速廓线图;(c)图(b)中黑色实线的平均Wigner-Ville分布
Fig. 1. Simulation results of tail vortex and Wigner-Ville distribution of radial velocity profiles. (a) Numerical simulation diagram of contour plot of tail vortex pair; (b) three radial velocity profiles of line of slight; (c) average Wigner-Ville distribution of black solid lines in fig. (b)
2014年,法国航空航天研究院报道了1.5 μm全光纤单频激光雷达的风速随距离变化的谱图结果[35]。谱图分析结果如
图 2. 风速随距离变化的谱图图像。(a)1.5 μm全光纤单频激光雷达;(b)长距多普勒测风雷达
Fig. 2. Spectral images of wind speed varying with distance. (a) 1.5 μm all-fiber single frequency lidar; (b) long-distance Doppler lidar
2016年,厦门大学提出采用EMD和极限学习机方法来准确预测风速序列[37]。该方法采用的风速序列是非线性非平稳序列,采用EMD将其分解为不同频率分量以减少非平稳性。为了避免极限学习机输入维数的随机性,重构各分量的相空间,建立各分量的极限学习模型预测风速序列,并对其进行叠加获得风速预测结果。与小波分解的风速预测结果对比表明,EMD和极限学习机的组合预测方法具有更高的判断精度。
2017年,美国科罗拉多大学Boulder分校提出一种小波自动识别方法,对美国南极麦克默多科考站的中间层及低热层的激光雷达观测数据进行分析,并对中间层和低热层重力波演化展开分析讨论[38]。小波自动识别方法提取出二维单色小波包,反演小波周期、波长、相速度、时间间隔和空间间隔等特性,其分析结果与2016年一维小波的分析结果一致[39]。
图 3. 二维小波重建结果。(a) 2014年7月16至18日相对温度扰动;(b)重建周期为3.6 h;(c)重建周期为4.8 h;(d)重建周期为7.8 h;(e)结合三个主小波重建的温度扰动场
Fig. 3. Reconstruction results of 2D wavelet. (a) Original relative temperature perturbations from July 16 to 18, 2014; (b) reconstruction period of 3.6 h; (c) reconstruction period of 4.8 h; (d) reconstruction period of 7.8 h; (e) the temperature perturbation field reconstructed from combining the above three major wave packets
图 4. 重力波扰动(a)~(c)和小波谱能量分布函数(d)~(f)。(a)原始温度扰动;(b)相位上行的波;(c)相位下行的波;(d)垂直波长与相速度;(e)垂直波长与周期;(f)高度与垂直波长
Fig. 4. Gravity wave perturbations (a)-(c) and distribution function of spectral energy (d)-(f). (a) Initial temperature perturbations; (b) waves with upward phase progression; (c) waves with downward phase progression; (d) Vertical wavelength versus phase velocity; (e) vertical wavelength versus period; (f) altitude versus vertical wavelength
2017年,法国航空航天研究院和法国放射性废物管理机构采用谱图方法,分析了1645 nm全光纤相干激光雷达在石油和天然气设施处对甲烷通量的测量能力[40]。
2017年,德国宇航局采用二维小波分析对南半球中纬度地区新西兰(45°S)和北半球高纬度地区芬兰(67°N)的两台瑞利激光雷达观测数据进行对比研究,观测30 km至85 km海拔高度冬季重力波传播特征[41]。
2018年,中国科学技术大学将自适应最优核时频分布应用于1.5 μm相干多普勒激光雷达实测数据,与传统的快速傅里叶方法对比风速反演结果[42]。
图 5. 模拟风切变与实测结果对比。(a)模拟结果;(b)实测
Fig. 5. Comparison of wind shear distribution between simulation results and actual measurements. (a) Simulation results; (b) actual measurements
2.2 信号去噪
噪声中未知信号的检测和提取是激光雷达信号处理中的难点。传统的方法是进行频率滤波或时间门选通,但是频率滤波不能将频率通带内的噪声滤除,时间门选通也不能把门限以内的噪声去除。通过时频变换,随机噪声的能量趋向于分布在整个时频域,而信号能量往往聚集在有限区间内,淹没在噪声中的信号在联合时频域更容易被检测出来。通过时变频率滤波还可以提高信噪比。
2003年,尼日利亚拉各斯大学和日本福冈大学采用短时傅里叶变换对平流层激光雷达后向散射信号进行时频分析,通过中值滤波器提取出淹没在强噪声中的瞬时脉冲,提高了激光雷达的信噪比[43]。
2006年,中国海洋大学采用EMD来提高Cabannes-Mie多普勒激光雷达信号的信噪比[44]。
图 6. 反演结果对比图。(a)原始数据和去噪信号;(b)去噪数据和1000组累加信号的平均值
Fig. 6. Comparison diagrams of inversion results. (a) Original and denoised data; (b) denoised data and average of 1000 sets of accumulative signals
2009年,哈尔滨工程大学基于重排小波-Radon变换,提高了线性调频雷达回波信号的聚集性,信噪比与未重排算法相比有3 dB左右的提升,避免了交叉项形成伪尖峰,削弱了噪声影响[45]。对于单分量信号以及信号幅度相近的多分量信号,该算法首先对线性调频信号进行小波变换,得到时间-尺度二维图像。将时间-尺度图按照映射关系转换为时间-频率图,对时间-频率图重排后进行Radon变换,设定阈值得到峰值点坐标,求取参数。对于强弱信号的多分量参数估计方法,采用上述步骤先估计出最强分量参数,将次强分量从信号中去除,然后估计次强分量,将其去除,依次类推。如果不确定信号数目,则根据剩余信号能量来决定算法是否继续。
2010年,西安理工大学采用EMD对激光雷达回波信号进行去噪处理[46]。EMD将激光雷达回波信号分解为不同尺度的本征模函数,而噪声往往聚集在高频分量中。采用Savitzky-Golay滤波器滤除高频分量,实验验证结果表明在保持有用信号下能够有效去除噪声。2011年,海**程大学根据激光雷达后向散射信号特性进行仿真生成,对固有模态分解滤波、扩展卡尔曼滤波和小波阈值滤波算法的滤波效果进行对比[47]。研究表明扩展卡尔曼滤波在均匀大气高信噪比时滤波效果较好;小波阈值滤波算法自适应性更强,在信噪比相同时表现最优;而固有模态分解滤波在高信噪比大气条件下滤波效果最优,可作为现有滤波算法的必要补充。2012年,中国科学院安徽光学精密机械研究所采用Hilbert-Huang变换辅助分析云层参数及垂直能见度并验证去噪效果[48]。对不同大气天气条件下的实验数据和反演结果进行对比发现,Hilbert-Huang时频方法在保留有用信息特征、有效抑制各种平稳或非平稳噪声基础上,降低了云高漏报率,提高了激光识别能力,增加了能见度反演精度。
2015年,美国陆军航空及导弹研发中心和美国国家航空航天局基于小波变换,对在海拔7000 inch(2133 m)的激光雷达采集的AS3500B直升机的叶片涡流交互噪声的实测数据进行时频分析[49]。与传统的傅里叶积分方法相比,小波变换在叶片涡流交互噪声处理方面取得了更好的聚集性效果。
2015年,西班牙加泰罗尼亚理工大学在世界统一时间(UTC)2013年4月20日的23时到23.5时时段,对高度为500 m的激光雷达原始后向散射信号进行谱图分析[50]。
2017年,中国科学院长春光学精密机械与物理研究所和中国科学院大学设计了激光多普勒测振仪(LDV)侦测喉咙振动,采用谱图方法和Wiener滤波后的谱图方法对光电探测器接收到的原始信号进行时频处理[51]。
2018年,西安理工大学提出基于小波阈值去噪算法对白天太阳光背景光进行滤波处理,实现对拉曼激光雷达回波信号的真实信号与噪声信号的分离[52]。
2018年,南京大学针对激光雷达信号非线性非稳定性特征,基于软阈值和粗糙惩罚的经验模态分解(EMD-STRP)对信号进行去噪处理[53]。有限长度的激光雷达数据经EMD分解后可得相同数据长度的固有模态函数,以噪声为主的IMF进行自相关运算后得到的自相关函数表现出噪声自相关的统计特性,以信号为主的IMF其自相关函数表现出普通信号的自相关统计特性。与基于相关EMD部分重建方法、硬阈值EMD、小波分析结果相比,EMD-STRP信噪比提高至22.25 dB,最远探测范围为11 km。
图 8. 谱图分析结果对比。(a)原始LDV谱图与波形图;(b) Wiener滤波后的信号谱图与波形图;(c)纯净信号的谱图与波形图
Fig. 8. Comparison of the spectrogram results. (a) Spectrogram and oscillogram of an original LDV signal; (b) spectrogram and oscillogram of a Wiener filtered signal; (c) spectrogram and oscillogram of a clean signal
图 9. 2016-09-22T00:00—2016-09-23T00:00,水汽混合比时间高度计(THI)图去噪前后结果显示图。 (a)去噪前; (b)去噪后
Fig. 9. THI displays of water-vapor mixing ratio recorded from 2016-09-22T00:00 to 2016-09-23T00:00 before and after denosing. (a) Before denoising; (b) after denoising
2.3 运动目标成像与检测
对激光雷达回波采用联合时频变换(JTFT),可以得到一个二维或三维的时频矩阵,其中可将信号时频谱视为时频图像进行分析。当分辨单元远小于目标尺寸时,选择高分辨率低交叉项干扰的时频变换,可以实现目标的高分辨率和多参数探测。
2003年,瑞典林雪平大学采用修正的果蝇优化算法(FOA)、Wigner-Ville分布和谱图分布进行目标特征识别[54]。修正的FOA对多普勒估计的输入参数敏感,而谱图和平滑Wigner-Ville分布的分析性能依赖于信号类型、平滑方法、窗函数长度等。
2005年,美国NASA采用全光纤相干多普勒激光雷达来协调行星精确软着陆,并对距离激光雷达250 m的目标回波信号进行谱图分析[55]。
图 10. 距离250 m处的目标回波信号谱图分析。(a)静目标;(b)动目标
Fig. 10. Spectrograms of the received signals from the targets at 250 m. (a) Stationary target; (b) moving target
2006年,美国爱达荷大学引入二维空间小波分析,根据激光雷达实测的数据分析阔针叶混合林的位置、高度和直径[56]。常用的可变窗方法需要树冠直径等先验知识,另外树冠和树高的弱相关性也会影响算法的灵敏性。而二维空间小波分析不需要雷达数据提供结构维度等先验知识,在检测植被结构和生物量方面颇具潜力。对于树高,小波分析优于可变窗分析方法;对于树冠,小波分析偏差较小,为-7%,可变窗分析方法偏差偏大,为-15%。
2006年,英国雷丁大学采用Gabor小波变换从机载激光雷达数据中分离地面和非地面[57]。Gabor小波具有最小的时频窗,可以在时域或频域获取理想的位置信息,且Gabor函数与人眼生物作用相仿,Gabor小波选取恰当参数可以很好地进行图像分割和纹理识别。Gabor小波捕捉局部不连续的响应,对于平坦区域其响应值低;而对于有房屋和植被的居民区域其响应值较高。将Gabor小波分布划分为
图 11. Gabor小波变换测试结果。(a)原始数据1;(b)原始数据1分割结果;(c)原始数据2;(d)原始数据2分割结果
Fig. 11. Test results of Gabor wavelet transform. (a) Tile 1 original data; (b) Tile 1 segmented result; (c) Tile 2 original data; (d) Tile 2 segmented result
2009年,美国南佛罗里达大学利用佛罗里达州大西洋激光雷达海岸线数据(该数据由美国海洋和大气管理局提供),基于小波变换和S变换对该激光雷达空间序列进行了时频分析,多尺度获取数千米的海岸线地貌演变特征[58]。实验结果与传统的三角间距测量方法一致,结合陆基雷达,该方法可以表征三角地貌和海岸线动态模型的演变。
2011年,美国橡树岭国家实验室采用压缩感知和正交匹配追踪算法从激光雷达LAS数据集中分离出建筑物和植物[59]。与传统信号处理方法不同,压缩感知理论对激光雷达数据进行稀疏分解,通过求解最优
图 12. 匹配追踪法分离树木和建筑物结果对比。(a)树木;(b)建筑物;(c)采用11×11窗检测到的树木区域;(d)采用11×11窗测到的建筑物区域;(e)采用7×7窗探测到的树木区域;(f)采用7×7窗探测到的建筑物区域
Fig. 12. Comparison of segmented trees and buildings using matching pursuit method. (a) Trees; (b) buildings; (c) tree area detected by an 11×11 window; (d) building area detected by an 11×11 window; (e) tree area detected by a 7×7 window; (f) building area detected by a 7×7 window
2012年,空**程大学针对逆合成孔径激光(ISAL)雷达回波信号的特点[60],利用重排Wigner-Ville分布和Hough变换对光外差探测后的信号进行时频分析来估计目标的运动速度,构造有效的补偿因子,完成回波信号的精确运动补偿,并进一步采用Keystone变换完成对目标散射点的越距离单元徙动的校正。该成像算法有效解决了激光信号极高载频、极大带宽和极短波长带来的问题,通过与微波波段逆合成孔径雷达的比较,证明了逆合成孔径成像激光雷达可实现对运动目标更快速、更高分辨的成像。
图 13. 谱图结果。(a)当纵模间隔为10 GHz时目标风速的归一化谱图;(b)硬阈值处理后的风速谱图
Fig. 13. Spectrogram results. (a) Normalized spectrogram of the target speed versus time with tone spacing of 10 GHz; (b) velocity spectrogram after hard threshold processing
2013年,爱沙尼亚塔林国立理工大学为了从高光谱激光雷达的诱导荧光光谱中提取信号特征,采用离散小波变换进行多尺度细化分析,来检测和分类水中的油污染[61]。基于稀疏优化方法,自动进行特征分类。分类结果和敏感性与一定数量的油污染离散通道传感器数据进行对比,模拟结果表明该方法准确性高,对单组对象中的有机化合物可进一步分类。
2013年,澳大利亚机器人中心将Velodyne激光雷达传感器置于移动车辆上收集数据,采用匹配追踪算法进行无监督场景识别分类[62],首次将深度特征学习应用于户外三维数据。所有目标采用的背景图像训练集是从Velodyne激光雷达传感器中随机选择的共20000组城市场景图像。分层匹配追踪算法提取的时频信息包含多尺度、多位置的特征信息,包括相位、相对方位、遮挡和姿态等。
2014年,北京装备学院针对ISAL雷达回波信号中存在色散效应和方位向多普勒时变的问题,采用分数阶傅里叶变换来补偿距离像色散中的谱峰分裂和展宽。在运动补偿后,结合CLEAN技术实现对机动目标的方位成像。仿真实验证明了该方法能够实现对运动目标的高分辨率、高聚焦性的ISAL成像[63]。
2015年,意大利比萨圣安娜高等大学基于归一化谱图方法,对双频激光雷达接收的后向散射信号进行了时频分析[64]。
2017年,美国路易斯安那州立大学在低流量条件下采用激光雷达和水下多波束测深仪对外滩地形进行详细测量,基于Hilbert-Huang变换对陆上地形变化粗糙度进行多尺度评估[65]。
2018年,中国科学院电子学研究所基于平滑伪Wigner-Ville分布,利用线性调频信号带宽为6 GHz、载波为1550 nm、脉冲重复频率为16.7 kHz的ISAL雷达成功实现了对户外1 km处的大小为20 cm×60 cm的飞机模型成像[66]。
图 14. 飞机模型以及基于两种方法的对比图。(a)石制飞机模型的光学图;(b) FFT(快速傅里叶变换)方法成像结果图;(c) FFT方法方位多视图;(d) JTFT方法方位多视图
Fig. 14. Airplane model and imaging results based on two methods. (a) Optical photo of the airplane model made of stone; (b) image result based on the FFT(fast Fourier transformation) method; (c) azimuth multilook result based on the FFT method; (d) azimuth multilook result based on the JTFT method
2.4 微多普勒特征分析
微多普勒现象展示了振动和旋转部位与目标交互作用的特性。运动目标如行人手臂和四肢摆动、人体心跳和胸腔振动、直升机的旋转桨、船只上的旋转天线、船只颠簸、电动机、装甲车、汽车车轮等,都会产生微多普勒效应。微动目标的微多普勒特征可以反映出目标精细特征,对于空间、地面目标、海中目标的探测识别具有重要意义。作为已有时域或频域方法的补充,微多普勒的联合时频表示可以提供振动或者旋转频率的时间信息,可以用来辨识感兴趣的目标。
2003年,美国帕森斯公司采用实验比较多种二次型时频方法对激光雷达回波信号的分析效果[67]。实验采用蒙特卡罗方法对振动旋转目标进行仿真,考虑由于振动、旋转和翻滚引起的散斑效应等因素,生成线性激光雷达回波信号和正弦调频回波信号,
表 1. 连续光相干激光雷达LO噪声峰值性能
Table 1. Approximate peak to LO noise performances for continuous wave coherent lidar
|
2007年,烟台大学采用短时傅里叶变换对目标振动引起的微多普勒效应进行了仿真,结果表明短时傅里叶变换能够从激光雷达回波信号中提取目标的微多普勒特征[68]。
2008年,法国格勒诺布尔第一大学和法国航空航天宇航局根据多普勒效应,采用敏感测速仪传感器和相干激光雷达对大楼振动频率进行谱图分析[69]。该工作在地震和结构工程力学方面具有广泛的应用前景。
2011年,空**程大学和复旦大学针对ISAL激光雷达回波信号特性,基于二值数学形态学腐蚀膨胀运算和推广Hough变换的方法有效地提取了目标微多普勒特征,并证明了ISAL激光雷达对cm或mm量级微动观测的有效性[70]。先对图像进行平滑处理,将其转化为二值图像,利用腐蚀膨胀方法提取微多普勒图像的边缘信息,这样每一条变细的正弦曲线的主要特征如频率、振幅和初相都与原曲线相同,而基线位置略有差别。对检测到的两正弦曲线基线加权平均可以得到更准确的基线位置,在此基础上采用推广Hough变换可将时频平面的曲线检测问题转化为参数空间(
2013年,空**程大学提出一种基于重排Wigner-Ville分布和波形熵的鸟类目标运动状态判别方法[71]。该算法首先基于压缩感知理论对鸟类目标ISAL雷达回波信号进行降采样,对降采样信号进行重排Wigner-Ville变换。通过波形熵函数检验每时刻所对应的重排Wigner-Ville图像中波形聚集效果来判别鸟类目标的运动状态。当波形熵结果小于给定门限时,鸟类目标处于滑翔状态;当波形熵结果大于给定门限时,鸟类目标处于振翅状态。分别处理鸟类目标的两种不同状态,获得滑翔状态下的鸟类目标高分辨二维图像和翅膀振动状态下带来的微多普勒曲线,完成鸟类目标的识别。该研究工作为鸟类目标的准确识别带来新的思路和技术支撑。
2015年,美国陆军研究实验室采用谱图方法对激光雷达数据成像,在微多普勒图像上的一系列模糊的人类动作中识别出手臂、腿和躯干等[72]。
2017年,中国人民解放军电子工程学院在分析直升机、螺旋桨飞机和喷气式飞机的微多普勒时频特征基础上,将平滑Wigner-Ville分布得到的时频图像生成灰度图像,利用图像处理方法进行去噪,通过分析时频图像中曲线差异,提取出时频图灰度共生矩阵(GLCM)和Tamura纹理特征,采用支持向量机算法对飞机目标进行进一步特征分类[73]。实验结果进一步验证了时频图纹理特征可以达到较为理想的飞机分类效果,利用改进的GLCM特征能够实现低信噪比条件下的目标准确分类。
3 时频分析在激光雷达中的应用趋势及挑战
时频分析方法是研究时变非平稳信号的有力工具,作为时间和频率的二维函数,时频分布给出了特定时间和特定频率范围的能量分布,也描述了非平稳信号的频率随时间变化的过程。当激光雷达信号呈现时变特性时,在联合时频域表示强度或能量分布的变换是最理想的变换。
时频分析在激光雷达大气参数提取、信号去噪、目标成像与检测、微多普勒分析等领域已有广泛的应用。时频分析在激光雷达中的应用有以下发展趋势。
1) 激光雷达的性能与分辨率密切相关。分辨率是指在描述目标向量时,激光雷达至少在某一维状态变量上对两个紧邻目标的区分能力。与目标相关的状态变量构成了一个四维空间,其中有方位角、俯仰角、距离和径向速度(径向速度即多普勒频移)。选择高分辨率低交叉项干扰的时频变换,可以实现目标的高分辨率和多参数探测。对大气参数的高分辨率多参数探测,有助于分析极端气象灾害,如大气湍流、飞机尾流和风切变等,在大气研究和****等领域具有重要意义。
2) 研究时频空间中的目标、杂波、传输通道中的区别,将在杂波抑制、传播效应消除和目标特征增强等方面发挥重要作用。比如海洋表面杂波干扰很大,这使漂浮小目标的检测和分类变得十分困难;在植被穿透应用中,穿透树冠的双向传播过程是隐蔽目标被检测到的重要因素。
3) 实际情况下,目标常常存在复杂的机动,包括平移和旋转运动,因此,要得到目标的聚焦图像就必须进行运动补偿。对于运动补偿中的多普勒跟踪、距离漂移和相位估计而言,时频分析是一种极具吸引力的方法。
4) 对于目标辨识应用,利用提取的时频特征提高分类器性能,有助于加强对复杂形状和奇异材料的散射机制的认识。没有时频分析,很难从测量数据中直观理解散射机制,并且目标的复杂性使解析求解非常困难。
5) 作为已有时域或频域方法的补充,时频分析在捕获时间相关的频率特征方面十分有用,可以给出微多普勒的时域和频域特征,因此,可以基于微多普勒现象识别感兴趣的目标。
6) 时频算法常与机器学习算法相结合,来提高激光雷达目标智能检测和认知能力。目前基于时频分析的感知算法,可实现对目标的分割、检测、跟踪、识别和重建等。
目前仍需要深入研究和优化时频算法,时频方法面临诸多挑战。
表 2. 各种时频方法对比
Table 2. Comparison of various time-frequency analysis methods
|
1) 与成像的傅里叶变换相比,基于时频的成像方法能够改善信噪比。但利用时频分析到底能改善多少信噪比仍属于未解问题,仍需对时频图像中信噪比进行进一步的定量分析。
2) 时频分析在运动补偿方面能够很好地与实际运动相结合,但是涉及更为复杂的运动时,如快速运动的空中目标和波涛起伏的海面上的船只,目标旋转不仅仅局限于二维平面,如何进行运动补偿是一个挑战性难题。
3) 另外一个问题是在激光雷达信号处理中时频方法的计算速度,例如自适应时频方法的缺点是当运动参数穷尽搜索时,计算代价通常很高,目前需要对一般基函数的快速时频方法需要做更多的研究工作。
4) 在时频变换获得的图像分辨率方面,仍需要寻找高时频分辨率和低交叉项干扰的快速时频变换。
4 结束语
时频分析在激光雷达信号的检测、估计、识别、分类、特征提取、信号分离、信号综合、滤波以及干扰、噪声的抑制等方面具有广泛的研究前景。联合时频分析是一种有效的信号和图像分析方法,因此,加大联合时频分析的研究力度有助于推动我国的空间技术、空间科学、空间应用步入先进行列。
[2] Wu SH, Zhai XC, Liu BY, et al. Characterization of aircraft dynamic wake vortices and atmospheric turbulence by coherent doppler lidar[C]∥The 28th International Laser Radar Conference, June 25-30, 2017, Bucharest, Romania. Amsterdam: EDP Sciences, 2018, 176: 06001.
[4] 麻晓敏, 陶宗明, 张璐璐, 等. 侧向散射激光雷达探测白天近地面气溶胶探测技术[J]. 光学学报, 2018, 38(4): 0401005.
[5] 叶光豪, 邓愫愫, 徐文兵, 等. 机载激光雷达技术应用于沙丘变形监测的研究[J]. 激光与光电子学进展, 2018, 55(5): 052802.
[8] BoashashB. Time-frequency signal analysis and processing: a comprehensive reference[M]. New York: Academic Press, 2015.
[9] Gabor D. Theory of communication. Part 1: The analysis of information[J]. Journal of the Institution of Electrical Engineers-Part III: Radio and Communication Engineering, 1946, 93(26): 429-441.
[10] Potter RK, Kopp GA, Green HC. Visible speech[M]. New York: Van Nostrand and Company, 1947.
[14] CohenL. Time-frequency analysis[M]. New Jersey: Prentice Hall PTR, 1995.
[15] BertrandJ, BertrandP. Affinetime-frequency distributions[C]∥1990 Special Conference on Time-Frequency Signal Analysis/International Symp on Signal Processing and its Applications, 1990, Gold Coast Australia. Melbourne: Longman Cheshire, 1992: 118- 140.
[19] Ville J. Theorie et application de la notion de signal analytique[J]. Cables et Transmission, 1948, 2(1): 61-74.
[21] Claasen T. Mecklenbrauker W F G. The Wigner distribution—a tool for time-frequency signal analysis[J]. Philips Journal of Research, 1980, 35(3): 217-250.
[22] FlandrinP, MartinW. A general class of estimators for the Wigner-Ville spectrum of non-stationary processes[M] ∥Bensoussan A, Lions J L. Analysis and Optimization of Systems. Berlin, Heidelberg: Springer, 1984: 15- 23.
[23] Born M, Jordan P. Zur quantenmechanik[J]. Zeitschriftfür Physik, 1925, 34(1): 858-888.
[26] Page C H. Instantaneous power spectra[J]. Journal of Applied Physics, 1952, 23(1): 103-106.
[30] Huang NE. Hilbert-Huang transform and its applications[M]. 2nd ed. New Jersey: World Scientific, 2014.
[34] BrousmicheS. Simulation of coherent Doppler LIDAR signals and their analysis with the Cohen's class: application to algorithms design for wake vortex detection and characterization[D]. Belgium: UCL-Université Catholique deLouvain, 2010.
[35] RenardW, GoularD, VallaM, et al. Beyond 10 km range wind-speed measurement with a 1.5 μm all-fiber laser source[C]∥2014 Conference on Lasers and Electro-Optics (CLEO)-Laser Science to Photonic Applications, June 8-13, 2014, San Jose, CA, USA. New York: IEEE, 2014: 14822367.
[37] Qiu JH, Shen SP, Xu GY. Short-term wind speed forecasting by combination of masking signal-based empirical mode decomposition and extreme learning machine[C]∥2016 11th International Conference on Computer Science & Education (ICCSE), August 23-25, 2016, Nagoya, Japan. New York: IEEE, 2016: 581- 586.
[45] 李利, 司锡才, 柴娟芳, 等. 基于重排小波-Radon变换的LFM雷达信号参数估计[J]. 系统工程与电子技术, 2009, 31(1): 74-77.
[46] Zhang YK, Ma XC, Hua DX, et al. An EMD-based denoising method for lidar signal[C]∥2010 3rd International Congress on Image and Signal Processing, October 16-18, 2010, Yantai, China. New York: IEEE, 2010: 4016- 4019.
[47] 陈冬, 王江安, 康圣. 脉冲激光雷达信号降噪方法对比[J]. 舰船科学技术, 2011, 33(4): 93-97.
[48] 何俊峰, 刘文清, 张玉钧, 等. HHT在激光云高仪后向散射信号处理中的应用[J]. 红外与激光工程, 2012, 41(2): 397-403.
[49] Stephenson JH, GreenwoodE. Effects of vehicle weight and true versus indicated airspeed on BVI noise during steady descending flight[C]∥71st Annual AHS Forum and Technology Display, May 5-7, 2015, Virginia Beach, VA, USA.2015.
[50] SaeedU, RocadenboschF, CrewellS. Adaptive estimation of the stable boundary-layer height using backscatter LiDAR data and a Kalman filter[C]∥2015 IEEE International Geoscience and Remote Sensing Symposium (IGARSS), July 26-31, 2015, Milan, Italy. New York: IEEE, 2015: 3591- 3594.
[51] Zhang H Y, Lv T, Yan C H. The novel role of arctangent phase algorithm and voice enhancement techniques in laser hearing[J]. Applied Acoustics, 2017, 126: 136-142.
[52] 王玉峰, 曹小明, 张晶, 等. 基于小波去噪算法的全天时大气水汽拉曼激光雷达探测与分析[J]. 光学学报, 2018, 38(2): 0201001.
[53] Chang J H, Zhu L Y, Li H X, et al. Noise reduction in Lidar signal using correlation-based EMD combined with soft thresholding and roughness penalty[J]. Optics Communications, 2018, 407: 290-295.
[54] Olsson, A. Target recognition by vibrometry with a coherent laser radar: LITH-ISY-EX-3050-2003[R/OL].( 2003-05-13)[2018-05-01]. http:∥www.ep.liu.se/exjobb/isy/2003/3050/.
[55] AmzajerdianF, PierrottetD, Tolson RH, et al. Development of a coherent LiDAR for aiding precision soft landing on planetary bodies[C]∥13th Coherent Laser Radar Conference, October 16-21, 2005, Kamakura, Japan.2005: 20050240846.
[56] Falkowski M J. Smith A M S, Hudak A T, et al. Automated estimation of individual conifer tree height and crown diameter via two-dimensional spatial wavelet analysis of lidar data[J]. Canadian Journal of Remote Sensing, 2006, 32(2): 153-161.
[57] WeiH, BartelsM. Unsupervised segmentation using Gabor wavelets and statistical features in LIDAR data analysis[C]∥18th International Conference on Pattern Recognition, August 20-24, 2006, Hong Kong, China. New York: IEEE, 2006: 667- 670.
[60] 何劲, 张群, 杨小优, 等. 逆合成孔径成像激光雷达成像算法[J]. 红外与激光工程, 2012, 41(4): 1094-1100.
[61] Sobolev I, Babichenko S. Analysis of the performances of hyperspectral lidar for water pollution diagnostics[J]. EARSeL eProceedings, 2013, 12(2): 113-123.
[62] Deuge MD, QuadrosA, HungC, et al. Unsupervised feature learning for classification of outdoor 3D scans[C]∥Proceedings of the2013Australasian Conference on Robitics and Automation, December 2-4, 2013, University of New South Wales, Sydney Australia. Australian Robotics and Automation Association, 2013, N/A:98586.
[63] Wu YH, RuanH, Yu DB. Inverse synthetic aperture laser radar imaging algorithm for maneuvering targets[C]∥2014 7th International Congress on Image and Signal Processing (CISP), October 14-16, 2014, Dalian, China. New York: IEEE, 2014: 569- 574.
[68] 王学勤, 董艳群, 原帅, 等. 激光雷达微多普勒效应的仿真研究[J]. 激光技术, 2007, 31(2): 117-119,146.
[70] 何劲, 张群, 罗迎, 等. 逆合成孔径成像激光雷达微多普勒效应分析及特征提取[J]. 电子学报, 2011, 39(9): 2052-2059.
[71] 朱丰, 张群, 冯有前, 等. 逆合成孔径激光雷达鸟类目标压缩感知识别方法[J]. 红外与激光工程, 2013, 42(1): 256-261.
[73] 王云鹏, 胡以华, 雷武虎, 等. 基于激光回波时频图纹理特征的飞机目标分类方法[J]. 光学学报, 2017, 37(11): 1128004.
Article Outline
刘燕平, 王冲, 夏海云. 时频分析在激光雷达中的应用进展[J]. 激光与光电子学进展, 2018, 55(12): 120005. Yanping Liu, Chong Wang, Haiyun Xia. Application Progress of Time-Frequency Analysis for Lidar[J]. Laser & Optoelectronics Progress, 2018, 55(12): 120005.