留言板

尊敬的读者、作者、审稿人, 关于本刊的投稿、审稿、编辑和出版的任何问题, 您可以本页添加留言。我们将尽快给您答复。谢谢您的支持!

姓名
邮箱
手机号码
标题
留言内容
验证码

星载激光测高仪多模式回波参数提取方法(特邀)

朱天豪 周辉 石岩 张千胤

朱天豪, 周辉, 石岩, 张千胤. 星载激光测高仪多模式回波参数提取方法(特邀)[J]. 红外与激光工程, 2022, 51(1): 20210836. doi: 10.3788/IRLA20210836
引用本文: 朱天豪, 周辉, 石岩, 张千胤. 星载激光测高仪多模式回波参数提取方法(特邀)[J]. 红外与激光工程, 2022, 51(1): 20210836. doi: 10.3788/IRLA20210836
Zhu Tianhao, Zhou Hui, Shi Yan, Zhang Qianyin. Parameter extraction method on the multiple mode waveforms of satellite laser altimeter(Invited)[J]. Infrared and Laser Engineering, 2022, 51(1): 20210836. doi: 10.3788/IRLA20210836
Citation: Zhu Tianhao, Zhou Hui, Shi Yan, Zhang Qianyin. Parameter extraction method on the multiple mode waveforms of satellite laser altimeter(Invited)[J]. Infrared and Laser Engineering, 2022, 51(1): 20210836. doi: 10.3788/IRLA20210836

星载激光测高仪多模式回波参数提取方法(特邀)

doi: 10.3788/IRLA20210836
基金项目: 国家自然科学基金(41971302);国防科工局十三五民用航天技术预先研究项目(250000887)
详细信息
    作者简介:

    朱天豪,男,硕士生,主要从事激光遥感与光电检测方面的研究

    周辉,男,教授,博士,主要从事激光遥感与空间探测等方面的研究

  • 中图分类号: TN958.98

Parameter extraction method on the multiple mode waveforms of satellite laser altimeter(Invited)

Funds: National Natural Science Foundation of China (41971302);The 13th Five-Year Civil Aerospace Technology Advance Research Project(250000887)
  • 摘要: 全波形星载激光测高仪的接收波形特征参数可以用于反演目标的形貌信息,传统的波形处理算法不能用于混叠严重以及偏离高斯形态的多模式波形特征参数提取。针对混叠严重的多模式回波,提出一种基于偏正态拟合模型,使用激励Richardson-Lucy反卷积算法、逐层分解算法、梯度下降法和非线性最小二乘拟合算法相组合的波形特征参数提取方法。采用已知参数的波形数据集、机载仿真波形数据集和全球生态系统动态调查(GEDI)激光雷达波形数据,基于波形相关系数与均方根误差(RMSE)、波形特征参数相对误差、波形分量个数提取正确率等评价指标开展波形处理试验,并将处理结果与传统的高斯分解结果进行比较分析。已知参数波形数据集处理结果的平均波形相关系数提升了约2%,RMSE降低了约47%,波形特征参数相对误差平均降低了约5%,分量个数提取正确率提升了约34%;机载仿真数据和GEDI波形数据处理结果的平均波形相关系数分别提升了约1%和2%,RMSE分别降低了约56%和54%。同时,开展了陡坡区域植被高度解算的仿真试验,得到的植被高度准确程度明显高于传统方法。所有处理结果均表明该方法更有利于多模式回波特征参数的提取以及目标参数的反演。
  • 图  1  总体算法流程图

    Figure  1.  Flow chart of the proposed method

    图  2  基于激励反卷积法提取波形分量峰值与个数的示意图

    Figure  2.  Sketch map of extracting the peak positions and numbers of the target response waveform (TRW) based on the boosted RL deconvolution method

    图  3  实验区域以及GEDI和机载激光雷达数据轨迹

    Figure  3.  Location of the study site and data trajectories of the GEDI and airborne lidar

    图  4  已知参数目标响应波形的特征参数及波形示例

    Figure  4.  Prescribed parameter distributions of the known-parameter TRW and an example of the TRW

    图  5  DRET方法和DGDL方法对接收波形的分解结果

    Figure  5.  Decomposed results of the received waveform by the DRET and DGDL methods

    图  6  DRET和DGDL方法对已知参数波形数据集的处理结果

    Figure  6.  The processed results of the known-parameter waveform data set obtained by the DRET and DGDL methods

    图  7  DRET和DGDL方法的波形特征参数提取误差和分量个数提取正确率

    Figure  7.  Extracted errors of the waveform parameters and the successful detection rates of the component numbers by the DRET and DGDL methods

    图  8  DRET方法和DGDL方法对机载仿真数据的分解结果

    Figure  8.  Decomposition results of airborne simulation data by the DRET and DGDL methods

    图  9  DRET和DGDL方法对GEDI数据的分解结果

    Figure  9.  Decomposition results of the GEDI data by the DRET and DGDL methods

    图  10  植被目标的点云分布

    Figure  10.  Point cloud distributions of the vegetated target

    图  11  不同树高下植被高度解算的误差分布

    Figure  11.  Distributions of the canopy height errors under different tree heights

    表  1  已知参数波形数据集的仿真参数

    Table  1.   Parameter setting of known-parameter waveform data set

    SourceAmplitude/VPeak position/nsPulse width/ns
    System response12015.6
    Target response0.2-1300-4005-15
    下载: 导出CSV

    表  2  两种方法提取结果的详细波形参数以及评价参数

    Table  2.   Detailed waveform parameters and evaluation indicators of the results extracted by the two methods

    MethodComponentIubaEvaluation criterion
    $ {C}_{x} $$ {\delta }_{x} $$ {\tau }_{I} $$ {\tau }_{u} $$ {\tau }_{b} $
    DRET10.69345.057.580.010.991.371.42%0.05%3.69%
    20.51366.2610.430.052.00%0.74%2.52%
    DGDL10.61343.897.16NA 10.972.6312.86%2.11%8.34%
    20.22361.2115.02NA 1NA 2NA 2NA 2
    30.31368.5510.98NA 118.97%2.31%6.22%
    1此处NA表示DGDL方法中不存在偏度系数。
    2此处NA表示该分量为错误分量,无法计算对应的参数误差。
    下载: 导出CSV

    表  3  两种方法的评估指标平均值

    Table  3.   Average value of the evaluation indexes of the two methods

    Evaluation indexDRETDGDL
    $ {C}_{x} $ 0.987 0.966
    $ {\delta }_{x} $ 1.217 2.305
    $ {\tau }_{I} $ 2.18% 7.56%
    $ {\tau }_{u} $ 0.52% 3.69%
    $ {\tau }_{b} $ 2.33% 6.47%
    $ {r}_{s} $ 98.70% 64.25%
    下载: 导出CSV

    表  4  DRET方法和DGDL方法对于机载仿真数据处理结果的相关系数($C_x $)和RMSE($\delta_x $)

    Table  4.   Correlation coefficient ($C_x $) and RMSE ($\delta_x $) of the DRET and the DGDL methods for the processing results of the airborne simulation data

    Evaluation indexDRETDGDL
    $ {C}_{x} $Max0.9990.991
    Min0.9370.872
    Average0.9950.986
    Standard deviation0.0120.034
    $ {\delta }_{x} $Max9.32614.273
    Min0.6071.988
    Average1.8693.516
    Standard deviation1.2912.936
    下载: 导出CSV

    表  5  DRET和DGDL方法对于GEDI数据处理结果的相关系数($C_x $)和RMSE ($\delta_x $)

    Table  5.   Correlation coefficient ($C_x $) and RMSE ($\delta_x $) of the DRET and the DGDL methods for the processing results of the GEDI data

    Evaluation indexDRETDGDL
    $ {C}_{x} $Max0.9990.982
    Min0.9390.854
    Average0.9930.977
    Standard deviation0.0370.062
    $ {\delta }_{x} $Max9.41215.504
    Min0.4493.075
    Average1.9534.248
    Standard deviation1.3633.712
    下载: 导出CSV
  • [1] Wagner W, Ullrich A, Ducic V, et al. Gaussian decomposition and calibration of a novel small-footprint full-waveform digitising airborne laser scanner [J]. ISPRS Journal of Photogrammetry and Remote Sensing, 2006, 60(2): 100-112. doi:  10.1016/j.isprsjprs.2005.12.001
    [2] Mallet C, Bretar F. Full-waveform topographic lidar: State-of-the-art [J]. ISPRS Journal of Photogrammetry and Remote Sensing, 2009, 64(1): 1-16. doi:  10.1016/j.isprsjprs.2008.09.007
    [3] Shi J, Menenti M, Lindenbergh R. Parameterization of surface roughness based on ICESat/GLAS full waveforms: A case study on the Tibetan Plateau [J]. Journal of Hydrometeorology, 2013, 14(4): 1278-1292. doi:  10.1175/JHM-D-12-0130.1
    [4] Bye I J, North P R J, Los S O, et al. Estimating forest canopy parameters from satellite waveform LiDAR by inversion of the FLIGHT three-dimensional radiative transfer model [J]. Remote Sensing of Environment, 2017, 188: 177-189. doi:  10.1016/j.rse.2016.10.048
    [5] Yang Chenchen, Xie Junfeng, Han Baomin, et al. Correlation analysis between ICESat/GLAS altimetry accuracy and echo waveform [J]. Applied Laser, 2020, 40(2): 9. (in Chinese)
    [6] Hermosilla T, Ruiz L A, Kazakova A N, et al. Estimation of forest structure and canopy fuel parameters from small-footprint full-waveform LiDAR data [J]. International Journal of Wildland Fire, 2013, 23(2): 224-233.
    [7] Wang X, Cheng X, Gong P, et al. Earth science applications of ICESat/GLAS: A review [J]. International Journal of Remote Sensing, 2011, 32(23): 8837-8864. doi:  10.1080/01431161.2010.547533
    [8] Mongus D, Žalik B. Parameter-free ground filtering of LiDAR data for automatic DTM generation [J]. ISPRS Journal of Photogrammetry and Remote Sensing, 2012, 67: 1-12. doi:  10.1016/j.isprsjprs.2011.10.002
    [9] Hofton M A, Minster J B, Blair J B. Decomposition of laser altimeter waveforms [J]. IEEE Transactions on Geoscience and Remote Sensing, 2000, 38(4): 1989-1996. doi:  10.1109/36.851780
    [10] Qin Y, Vu T T, Ban Y. Toward an optimal algorithm for LiDAR waveform decomposition [J]. IEEE Geoscience and Remote Sensing Letters, 2011, 9(3): 482-486.
    [11] Zhu J, Zhang Z, Hu X, et al. Analysis and application of LiDAR waveform data using a progressive waveform decomposition method [J]. International Archives of Photogrammetry, Remote Sensing and Spatial Information Sciences, 2011, W12(5): 31-36.
    [12] Jutzi B, Stilla U. Range determination with waveform recording laser systems using a Wiener filter [J]. ISPRS Journal of Photogrammetry & Remote Sensing, 2007, 61(2): 95-107.
    [13] Xie Junfeng, Yang Chenchen, Mei Yongkang, et al. Full waveform decomposition of spaceborne laser based on genetic algorithm [J]. Infrared and Laser Engineering, 2020, 49(11): 20200945. (in Chinese)
    [14] Miller S D, Stephens G L. Multiple scattering effects in the lidar pulse stretching problem [J]. Journal of Geophysical Research: Atmospheres, 1999, 104(D18): 22205-22219. doi:  10.1029/1999JD900481
    [15] Zhang Z, Xie H, Tong X, et al. A combined deconvolution and Gaussian decomposition approach for overlapped peak position extraction from large-footprint satellite laser altimeter waveforms [J]. IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing, 2020, PP(99): 1-1.
    [16] Zhou T, Popescu S C, Krause K, et al. Gold – A novel deconvolution algorithm with optimization for waveform LiDAR processing [J]. ISPRS Journal of Photogrammetry and Remote Sensing, 2017, 129(1): 131-150.
    [17] Miroslav Morháč, Vladislav Matoušek. High-resolution boosted deconvolution of spectroscopic data [J]. Journal of Computational & Applied Mathematics, 2011, 235(6): 1629-1640.
    [18] Zhao Quanhua, Chen Weiduo, Wang Yu, et al. Variable component waveform decomposition of partial normal full wave lidar data [J]. Optics and Precision Engineering, 2018, 26(1): 161-171. (in Chinese)
    [19] 李勇, 范承玉, 时东锋. 基于加速正则化RL算法的大气湍流退化图像盲复原方法[J]. 大气与环境光学学报, 2011, 6(5): 9.

    Li Yong, fan Chengyu, Shi Dongfeng. Blind restoration method of atmospheric turbulence degraded image based on accelerated regularization RL algorithm [J] Journal of Atmospheric and Environmental Optics, 2011, 6 (5): 342-350. (in Chinese)
    [20] Morháč M. Deconvolution methods and their applications in the analysis of γ-ray spectra [J]. Nuclear Instruments and Methods in Physics Research Section A: Accelerators, Spectrometers, Detectors and Associated Equipment, 2006, 559(1): 119-123.
    [21] Luo Min, Shi Yan, Zhou Hui, et al. Lidar pulse waveform decomposition based on variable component parameter random sampling [J]. Infrared and Laser Engineering, 2019, 48(10): 1005009. (in Chinese)
    [22] Popescu S C, Zhao K, Neuenschwander A, et al. Satellite lidar vs. small footprint airborne lidar: Comparing the accuracy of aboveground biomass estimates and forest structure metrics at footprint level [J]. Remote Sensing of Environment, 2011, 115(11): 2786-2797. doi:  10.1016/j.rse.2011.01.026
    [23] Hudnut K W, Brooks B A, Scharer K, et al. Airborne lidar and electro‐optical imagery along surface ruptures of the 2019 Ridgecrest earthquake sequence, Southern California [J]. Seismological Research Letters, 2020, 91(4): 2096-2107. doi:  10.1785/0220190338
    [24] Hancock S, Armston J, Hofton M, et al. The GEDI simulator: A large‐footprint waveform lidar simulator for calibration and validation of spaceborne missions [J]. Earth and Space Science, 2019, 6(2): 294-310. doi:  10.1029/2018EA000506
    [25] Liu Ren, Xie Junfeng, Mo Fan, et al. Simulation of echo waveform of spaceborne laser altimeter based on fine terrain [J]. Acta Photonica Sinica, 2018, 47(11): 1128004. (in Chinese)
  • [1] 李莹莹, 刘子维, 张静坤, 吴琳琳, 纪雪, 王明常.  基于波形形态特征的单频机载激光雷达测深全波形数据分类 . 红外与激光工程, 2023, 52(9): 20230096-1-20230096-12. doi: 10.3788/IRLA20230096
    [2] 吴鹏飞, 邓植中, 雷思琛, 谭振坤, 王姣.  基于激光散斑图像多特征参数的表面粗糙度建模研究 . 红外与激光工程, 2023, 52(12): 20230348-1-20230348-11. doi: 10.3788/IRLA20230348
    [3] 梅永康, 谢俊峰, 陈伟, 刘仁.  多特征参数约束的星载激光高程控制点提取 . 红外与激光工程, 2022, 51(9): 20210997-1-20210997-12. doi: 10.3788/IRLA20210997
    [4] 孙俊灵, 马鹏阁, 郭清源, 韩红印, 李伟, 陶然.  低信噪比下机载平台多脉冲激光测距机目标回波降噪算法 . 红外与激光工程, 2021, 50(6): 20210005-1-20210005-9. doi: 10.3788/IRLA20210005
    [5] 陈宏宇, 罗海波, 惠斌, 常铮.  采用多特征融合的子块自动提取方法 . 红外与激光工程, 2021, 50(8): 20200407-1-20200407-9. doi: 10.3788/IRLA20200407
    [6] 左志强, 唐新明, 李国元, 李松.  GF-7星载激光测高仪全波形自适应高斯滤波 . 红外与激光工程, 2020, 49(11): 20200251-1-20200251-11. doi: 10.3788/IRLA20200251
    [7] 郭金权, 李国元, 左志强, 张宁, 裴亮, 卢刚.  高分七号卫星激光测高仪全波形数据质量及特征分析 . 红外与激光工程, 2020, 49(S2): 20200387-20200387. doi: 10.3788/IRLA20200387
    [8] 刘向锋, 黄庚华, 张志杰, 王凤香, 舒嵘.  高分七号激光测高中全波形回波数据的EMD降噪 . 红外与激光工程, 2020, 49(11): 20200261-1-20200261-10. doi: 10.3788/IRLA20200261
    [9] 谷牧, 任栖锋, 廖胜, 周金梅, 赵汝进.  基于点目标特征参数提取的红外多光谱设计 . 红外与激光工程, 2020, 49(5): 20190462-20190462-10. doi: 10.3788/IRLA20190462
    [10] 赵丽娟, 徐志钮, 李永倩.  布里渊增益谱特征参数提取准确性影响因素分析 . 红外与激光工程, 2019, 48(2): 222003-0222003(8). doi: 10.3788/IRLA201948.0222003
    [11] 张文豪, 李松, 张智宇, 刘芮, 马跃.  利用波形匹配实现卫星激光测高脚点精确定位的方法 . 红外与激光工程, 2018, 47(11): 1117007-1117007(8). doi: 10.3788/IRLA201847.1117007
    [12] 徐志钮, 胡志伟, 赵丽娟, 杨志, 陈飞飞, 李永倩, 陈永辉.  采用Voigt模型的布里渊散射谱关键特征高精度提取方法 . 红外与激光工程, 2018, 47(S1): 74-81. doi: 10.3788/IRLA201746.S122004
    [13] 马跃, 张文豪, 张智宇, 马昕, 李松.  基于半解析模型的激光测高回波海水海冰波形分类方法 . 红外与激光工程, 2018, 47(5): 506005-0506005(7). doi: 10.3788/IRLA201847.0506005
    [14] 孙俊灵, 马鹏阁, 孙光民, 羊毅.  低信噪比下机载多脉冲激光雷达姿态不敏感性特征提取研究 . 红外与激光工程, 2017, 46(3): 330002-0330002(9). doi: 10.3788/IRLA201746.0330002
    [15] 李少辉, 周辉, 倪国强.  基于星载激光测高仪多模式回波的激光测距修正值分析 . 红外与激光工程, 2017, 46(10): 1006001-1006001(8). doi: 10.3788/IRLA201759.1006001
    [16] 赵丽娟, 李永倩, 徐志钮.  Lorentzian型布里渊频谱特征提取时模型的影响 . 红外与激光工程, 2016, 45(5): 522002-0522002(6). doi: 10.3788/IRLA201645.0522002
    [17] 孙俊灵, 马鹏阁, 孙光民, 金秋春, 羊毅.  基于目标波形模型的多脉冲激光雷达目标信号模拟 . 红外与激光工程, 2016, 45(7): 726006-0726006(7). doi: 10.3788/IRLA201645.0726006
    [18] 胡欣, 张文攀, 殷瑞光, 李慧, 赵琳锋, 刘艳芳.  激光在锥形多模光纤中的耦合效率与传输模式 . 红外与激光工程, 2013, 42(2): 372-375.
    [19] 纪超, 刘慧英, 邵刚, 孙景峰.  基于生物激励计算模型在图像显著性提取中的研 . 红外与激光工程, 2013, 42(3): 823-828.
    [20] 李姣, 张天序.  基于Laplacian-Markov先验数据的加权光谱反卷积模型 . 红外与激光工程, 2013, 42(12): 3464-3469.
  • 加载中
图(11) / 表(5)
计量
  • 文章访问数:  334
  • HTML全文浏览量:  71
  • PDF下载量:  44
  • 被引次数: 0
出版历程
  • 收稿日期:  2021-11-10
  • 修回日期:  2021-12-29
  • 网络出版日期:  2022-02-10
  • 刊出日期:  2022-01-31

星载激光测高仪多模式回波参数提取方法(特邀)

doi: 10.3788/IRLA20210836
    作者简介:

    朱天豪,男,硕士生,主要从事激光遥感与光电检测方面的研究

    周辉,男,教授,博士,主要从事激光遥感与空间探测等方面的研究

基金项目:  国家自然科学基金(41971302);国防科工局十三五民用航天技术预先研究项目(250000887)
  • 中图分类号: TN958.98

摘要: 全波形星载激光测高仪的接收波形特征参数可以用于反演目标的形貌信息,传统的波形处理算法不能用于混叠严重以及偏离高斯形态的多模式波形特征参数提取。针对混叠严重的多模式回波,提出一种基于偏正态拟合模型,使用激励Richardson-Lucy反卷积算法、逐层分解算法、梯度下降法和非线性最小二乘拟合算法相组合的波形特征参数提取方法。采用已知参数的波形数据集、机载仿真波形数据集和全球生态系统动态调查(GEDI)激光雷达波形数据,基于波形相关系数与均方根误差(RMSE)、波形特征参数相对误差、波形分量个数提取正确率等评价指标开展波形处理试验,并将处理结果与传统的高斯分解结果进行比较分析。已知参数波形数据集处理结果的平均波形相关系数提升了约2%,RMSE降低了约47%,波形特征参数相对误差平均降低了约5%,分量个数提取正确率提升了约34%;机载仿真数据和GEDI波形数据处理结果的平均波形相关系数分别提升了约1%和2%,RMSE分别降低了约56%和54%。同时,开展了陡坡区域植被高度解算的仿真试验,得到的植被高度准确程度明显高于传统方法。所有处理结果均表明该方法更有利于多模式回波特征参数的提取以及目标参数的反演。

English Abstract

    • 全波形星载激光测高仪是一种具备接收脉冲波形记录功能的主动式遥感设备,其接收脉冲信号是目标响应函数和星载激光测高仪系统响应函数的卷积,它与被测目标的几何和物理属性有关。在大尺度的激光足印范围内,复杂的目标分布将使得接收波形通常呈现多模式形态分布[1-2]。通过对接收脉冲波形特征参数的提取可以实现陆地的粗糙度与坡度、冰层高度、森林冠层结构和地物反射率等信息的反演[3]。因此,全波形星载激光测高仪在地物目标分类、高精度立体测绘、森林调查和全球冰川消融监测等方面具有广泛的应用[4]

      全波形数据特征参数的提取精度直接决定星载激光测高仪数据应用的准确性和可靠性[5]。目前针对全波形激光测高仪数据,最常用的波形分解方法是直接高斯分解法,该方法基于高斯模型和拐点法对波形的初始参数进行估计,基于波形拟合优化算法实现波形特征参数的提取[6-8]。常用的波形拟合优化算法包括非线性最小二乘拟合Levenberg-Marquardt(LM)算法、最大期望的极大似然法、可逆马尔可夫链蒙特卡洛算法和遗传算法[9-13]。前两种算法对初始参数的精度要求较高,可逆马尔可夫链蒙特卡洛和遗传算法的时间效率较差,波形特征参数的提取结果不稳定。对于具有一定混叠程度的多模式回波,这些算法往往会陷入局部最优解,无法准确甚至错误地分解出波形的特征参数。

      反卷积算法可以有效消除发射脉冲和接收系统响应对接收脉冲波形的展宽效应,只保留反映被测目标特性的目标响应波形,从而有利于混叠波形分量的分解。常用的反卷积算法分为两大类:直接反卷积算法和迭代反卷积算法。直接反卷积算法包括维纳滤波算法和Tikhonov正则化滤波算法,它们所得到的目标响应波形容易受到噪声的影响并且稳定性较差。迭代反卷积算法包括Richardson-Lucy(RL)算法、Gold算法和盲去卷积算法[14-17],通过设置合适的迭代参数可以得到相对准确的目标响应波形。然而,迭代反卷积算法仍然存在三个方面的问题:一是目标响应波形可能会偏离标准高斯形态而出现拖尾分布;二是目标响应波形无法完全与真实目标响应波形一致;三是目标响应波形仍然会存在混叠现象,这些问题都将影响目标响应波形参数的提取精度。

      针对混叠严重以及偏离高斯形态的星载激光测高仪接收波形,提出一种基于偏正态模型,采用激励反卷积算法、逐层分解算法、梯度下降法和LM算法相组合的波形分解方法(DRET)。利用已知波形参数的数据集、机载仿真波形数据集和全球生态系统动态调查(GEDI)波形数据,实现波形处理算法性能的综合评估。

    • 星载激光测高仪接收脉冲波形特征参数的提取算法主要分为三个过程:使用带有激励的反卷积算法对接收波形进行反卷积运算,得到目标响应波形;基于逐层分解、梯度下降法、LM算法和偏正态分解模型提取目标响应波形的初始参数;采用卷积运算和目标响应波形的初始参数,结合LM算法提取接收波形的特征参数。接收脉冲波形特征参数提取的总体流程如图1所示。

      图  1  总体算法流程图

      Figure 1.  Flow chart of the proposed method

    • 激光足印范围内的复杂地物通常会使得接收脉冲回波的分量出现拖尾分布,偏离标准的高斯模型。因此,采用偏正态模型来描述接收脉冲信号的理论分布,该模型既可以兼容高斯模型,又能体现接收波形的拖尾效应。偏正态模型的理论表达式为:

      $$ \begin{split} Y\left(t,\boldsymbol{P}\right)=&{\displaystyle\sum }_{j=1}^{n}f\left(t,{\boldsymbol{p}}_{j}\right)+N=\\ &{\displaystyle\sum }_{j=1}^{n}\left[2{I}_{j}\varphi \left(\frac{t-{u}_{j}}{{b}_{j}}\right) \varPhi \left({a}_{j}\frac{t-{u}_{j}}{{b}_{j}}\right)\right]+N\ \end{split} $$ (1)

      式中:$ \boldsymbol{P}=\left[{\boldsymbol{p}}_{1},{\boldsymbol{p}}_{2},\cdots ,{\boldsymbol{p}}_{n}\right] $表示所有的波形分量的特征参数集合;$ {\boldsymbol{p}}_{j}=\left[{I}_{j},{u}_{j},{b}_{j},{a}_{j}\right] $表示每一个波形分量的特征参数,$ {I}_{j} $$ {u}_{j} $$ {b}_{j} $$ {a}_{j} $分别表示第$ j $个偏正态分量的幅值、峰值时刻、脉冲宽度以及偏度;$ f\left(t,{\boldsymbol{p}}_{j}\right) $表示波形分量;$ N $表示叠加在波形上的噪声;$ n $为波形分量的个数;$ \varphi \left(t\right) $表示标准高斯函数;$ \varPhi \left(t\right) $表示标准正态分布的概率密度分布函数。当偏度系数为0时,偏正态模型就变成了高斯模型[18]

      (1)目标响应波形的解算

      基于贝叶斯理论的RL算法是一种时间域上的迭代反卷积算法[19],其迭代计算模型为:

      $$ {m}^{i+1}\left(t\right)={m}^{i}\left(t\right)\cdot \left(\frac{{R}_{x}\left(t\right)}{{m}^{i}\left(t\right)\mathrm{*}h\left(t\right)}\mathrm{*}h\left(-t\right)\right) $$ (2)

      式中:$ {m}^{i}\left(t\right) $表示经过$ i $次迭代后的目标响应波形;$ {R}_{x}\left(t\right) $$ h\left(t\right) $分别表示接收脉冲信号与系统响应;*表示卷积运算。

      RL反卷积算法输出的目标响应波形与迭代次数有关,当接收波形各分量之间的混叠程度进一步加重时,仅依靠增加RL反卷积算法的迭代次数无法分离混叠的波形分量。因此,文中对反卷积的输出$ {m}^{i}\left(t\right) $进行非线性激励,然后再次进行迭代反卷积运算,从而进一步消除混叠效应。研究表明非线性幂运算激励最有效,其表达式为:

      $$ {m}_{k}\left(t\right)={\left[{m}^{i}\left(t\right)\right]}^{k} $$ (3)

      式中:$ k $为激励系数,取值范围在1~2之间[20]。通过查找经过激励的RL反卷积运算的目标响应波形的峰值位置及个数,提取每个波形分量的峰值时刻初始值$ {u}_{j}\left(j=\mathrm{1,2},\cdots ,n\right) $$ n $为波形分量个数的初始值,其原理如图2所示。

      图  2  基于激励反卷积法提取波形分量峰值与个数的示意图

      Figure 2.  Sketch map of extracting the peak positions and numbers of the target response waveform (TRW) based on the boosted RL deconvolution method

      (2)目标响应波形的分解

      逐层分解方法将原始波形中的子脉冲按照顺序逐个估计并从原始波形中分离,从而达到完整分解波形的效果。当原始波形中的最大子脉冲分离之后,剩余波形中产生新的最大子脉冲并可以进一步估计和分解。

      将由激励的RL反卷积算法输出的目标响应波形记为$ {O}_{t} $,根据目标响应波形峰值时刻的初始值$ {u}_{j} $,获取每个波形分量对应的幅值$ {I}_{j} $。从每一个峰值时刻$ {u}_{j} $向左右两边寻找幅值为88.25%$ {I}_{j} $的位置$ {u}_{j}^{L} $$ {u}_{j}^{R} $,两者的差值作为脉宽宽度的初始值$ {b}_{j} $,所有波形分量的偏度系数$ {a}_{j} $暂时设置为0。

      使用梯度下降法对上述初始参数进行优化,设定优化的目标函数为:

      $$ {O}_{p}={\sum }_{t={u}_{_L}}^{{u}_{_R}}{\left[{O}_{t}-f\left(t,{\boldsymbol{p}}_{j}\right)\right]}^{2} $$ (4)

      式中:$ {O}_{p} $表示在$ {u}_{j}^{L} $$ {u}_{j}^{R} $范围内偏正态分量与目标响应波形$ {O}_{t} $的残差平方和。迭代更新后的波形分量参数可以表示为:

      $$ \left\{\begin{array}{l}\begin{array}{l}{I}_{j}^{GD}={I}_{j}-{\alpha }_{I}\dfrac{\partial {O}_{p}}{\partial {I}_{j}}\\ {u}_{j}^{GD}={u}_{j}\end{array}\\ {b}_{j}^{GD}={b}_{j}-{\alpha }_{b}\dfrac{\partial {O}_{p}}{\partial {b}_{j}}\\ {a}_{j}^{GD}={a}_{j}-{\alpha }_{a}\dfrac{\partial {O}_{p}}{\partial {a}_{j}}\end{array}\right. $$ (5)

      式中:$ {\alpha }_{I} $$ {\alpha }_{b} $$ {\alpha }_{a} $表示不同参数的迭代步长,可以由Armijo-Goldstein准则和Wolfe-Powell准则计算得到。如果更新后的波形分量与当前的波形分量之间残差的标准差小于4.5倍的噪声水平[21],则停止迭代,得到优化后的第$ j $个偏正态分量,其参数表示为$ {\boldsymbol{p}}_{j}^{GD}=\left[{I}_{j}^{GD},{u}_{j}^{GD},{b}_{j}^{GD},{a}_{j}^{GD}\right] $。否则继续迭代更新,直到迭代步长小于10-8

      将第$ j $个分量从目标响应波形$ {O}_{t} $中分离得到新的目标响应波形,即为:

      $$ {O}_{t}={O}_{t}-f\left(t,{\mathit{p}}_{j}^{GD}\right)$$ (6)

      当新的目标响应波形$ {O}_{t} $的最大值大于噪声水平的4.5倍时,重复上述步骤,以提取所有波形分量的参数。将这些参数集合作为目标响应波形的初始值,基于LM优化算法,对目标响应函数$ {m}^{i}\left(t\right) $进行进一步的处理,最终提取得到目标响应波形的参数结果$ {\boldsymbol{X}}_{{\rm{all}}}=\left[{\boldsymbol{p}}_{1}^{{\rm{LM}}},{\boldsymbol{p}}_{2}^{{\rm{LM}}},\cdots ,{\boldsymbol{p}}_{n}^{{\rm{LM}}}\right] $

      (3)接收波形特征参数的提取

      考虑到接收波形是目标响应波形和系统响应的卷积,则接收波形的初始参数可以由目标响应波形的参数结果来表示:

      $$ \left\{\begin{array}{l}{I}_{j}^{Rx}=\dfrac{\sqrt{2\pi }{\cdot I}_{j}^{{\rm{LM}}}\cdot {I}^{h}\cdot {b}_{j}^{{\rm{LM}}}\cdot {b}^{h}}{\sqrt{{{b}_{j}^{{\rm{LM}}}}^{2}+{{b}^{h}}^{2}}}\\ {u}_{j}^{Rx}={u}_{j}^{{\rm{LM}}}\\ \begin{array}{l}{b}_{j}^{Rx}=\sqrt{{{b}_{j}^{{\rm{LM}}}}^{2}+{{b}^{h}}^{2}}\\ {{a}_{j}^{Rx}=a}_{j}^{{\rm{LM}}}\end{array}\end{array}\right. $$ (7)

      式中:$ {I}_{j}^{{\rm{LM}}} $$ {u}_{j}^{{\rm{LM}}} $$ {b}_{j}^{{\rm{LM}}} $$ {a}_{j}^{{\rm{LM}}} $表示第$ j $个目标响应波形分量的特征参数;$ {I}^{h} $$ {b}^{h} $分别表示系统响应波形$ h\left(t\right) $的幅值和脉冲宽度。

      接收波形分量需要采取去除和合并步骤,以保证最终分解得到的波形分量符合实际的接收回波特征。首先去除幅值小于4.5倍噪声水平以及面积小于相邻分量面积的5%的波形分量,然后将峰值时刻间隔小于2 ns的相邻波形分量进行合并,再基于偏正态模型和LM优化算法对接收波形进行拟合,以提取得到接收波形中所有分量的特征参数。

    • 采用三套指标体系来评估接收波形的处理效果:波形相关系数和均方根误差(RMSE)、波形特征参数的提取误差和波形分量提取正确率。

      (1)波形相关系数和均方根误差

      $$ \begin{array}{c}{C}_{x}=\dfrac{\displaystyle\sum \left({R}_{a}-\stackrel-{{R}_{a}}\right)\left({R}_{b}-\stackrel-{{R}_{b}}\right)}{\sqrt{\displaystyle\sum {\left({R}_{a}-\stackrel-{{R}_{a}}\right)}^{2}\displaystyle\sum {\left({R}_{b}-\stackrel-{{R}_{b}}\right)}^{2}}}\\ {\delta }_{x}=\dfrac{1}{{\delta }_{{\rm{noise}}}}\sqrt{\dfrac{\displaystyle\sum {\left({R}_{a}-{R}_{b}\right)}^{2}}{N-1}}\end{array} $$ (8)

      式中:$ {R}_{a} $表示原始波形;$ {R}_{b} $表示拟合后的波形;$ \overline{{R}_{a}} $$ \overline{{R}_{b}} $分别代表它们的均值;$ {\delta }_{{\rm{noise}}} $表示噪声的标准差。

      (2)波形分量特征参数的提取误差

      将提取得到每个波形分量的特征参数与其设定的特征参数进行差值运算,得到特征参数的相对误差为:

      $$ {\tau }_{\mu }=\frac{\left|{\mu }_{ex}-{\mu }_{kn}\right|}{{\mu }_{kn}}\times 100\text{%} $$ (9)

      式中:$ {\mu }_{ex} $表示提取得到的波形特征参数;$ {\mu }_{kn} $表示设定的波形特征参数。

      (3)波形分量提取正确率

      接收波形处理方法所提取得到的波形分量个数可能与实际波形分量个数存在差异,出现波形分量漏判的可能,定义波形分量提取正确率为:

      $$ {r}_{s}=\dfrac{{N}_{s}}{{N}_{T}}\times 100\text{%} $$ (10)

      式中:$ {N}_{s} $表示提取得到的波形分量个数;$ {N}_{T} $表示设定的波形分量个数。

    • 机载仿真波形数据集和GEDI实测数据集的实验区域选在美国加州里奇克莱斯特境内,其经纬度范围为[35°35′50′N-35°38′55″N,117°36′40″W-117°40′30″W],面积约为30 km2,数据轨迹如图3所示。实验区域内有丰富的地物种类,包括矮小灌木、城市建筑、人造公园以及起伏的矮坡等,这使得接收波形的分布形态更为多样化,从而有利于波形处理算法性能的验证。

      图  3  实验区域以及GEDI和机载激光雷达数据轨迹

      Figure 3.  Location of the study site and data trajectories of the GEDI and airborne lidar

    • 已知参数的波形数据集是由GEDI的系统响应与所设定的目标响应波形的卷积结果得到,每个目标响应波形都由两个高斯分量组成,高斯分量的幅值、峰值时刻和脉冲宽度均在一定范围内随机选取,详细的参数设置如表1所示。文中设置了2000组已知参数的目标响应波形,其特征参数分布和波形示例如图4所示。噪声水平会对接收波形的分解产生一定的影响,参考GEDI接收波形的信噪比数据,将所有已知参数波形数据集的信噪比统一设置为15 dB。

      表 1  已知参数波形数据集的仿真参数

      Table 1.  Parameter setting of known-parameter waveform data set

      SourceAmplitude/VPeak position/nsPulse width/ns
      System response12015.6
      Target response0.2-1300-4005-15

      图  4  已知参数目标响应波形的特征参数及波形示例

      Figure 4.  Prescribed parameter distributions of the known-parameter TRW and an example of the TRW

    • 机载仿真波形数据集是利用高精度的机载激光雷达数据和波形模拟方法得到的[22]。所选区域的机载激光雷达数据来源于Optech Titan多光谱激光雷达系统在2019年所采集的离散点云[23]。该系统的光束发散角为0.5 mrad,发射激光脉冲的重复频率为300~500 kHz,激光足印点的水平和高程误差分别小于45 cm和15 cm。激光雷达的点云密度达到33个/m2,这能更加真实地反映仿真目标的空间分布。

    • GEDI激光测高仪系统包括三个激光器,每个激光器的重复频率为242 Hz,能够同时产生八条相互平行的激光足印轨迹[24]。相邻足印轨迹在沿轨方向和垂轨方向上分别相距60 m和600 m。所选实验区域的GEDI数据采集于2019年8月,包含1000个有效的激光足印。GEDI的接收波形以及激光足印的准确位置可以在GEDI 1 B级产品中获取,波形特征参数的提取结果存储在GEDI 2 A级产品中。

    • 基于拐点法与LM算法相结合的直接高斯分解法(DGDL)已成功应用于美国的ICESat-1和GEDI星载激光测高仪的波形数据处理。文中将DRET方法与DGDL方法所提取得到的波形结果进行比较,以验证处理算法的精度和可靠性。

    • 选取一组具有一定混叠程度的接收波形作为处理对象,采用DRET与DGDL两种方法进行处理,其结果如图5所示。

      图  5  DRET方法和DGDL方法对接收波形的分解结果

      Figure 5.  Decomposed results of the received waveform by the DRET and DGDL methods

      图5中可以看出,原始波形的下降沿呈现一定的拖尾分布,DRET方法能够将接收波形成功分解成两个分量,在波形的峰值位置以及尾部都有很好的拟合效果。然而DGDL方法提取得到一个额外的错误分量,使得拟合波形在上升沿、峰值位置和下降沿出现一定程度的误差。高斯模型无法准确拟合这种拖尾分布的波形,而偏正态模型具有较好的拟合效果。将两种方法所提取得到的波形特征参数以及定量评价指标列在表2中。

      表 2  两种方法提取结果的详细波形参数以及评价参数

      Table 2.  Detailed waveform parameters and evaluation indicators of the results extracted by the two methods

      MethodComponentIubaEvaluation criterion
      $ {C}_{x} $$ {\delta }_{x} $$ {\tau }_{I} $$ {\tau }_{u} $$ {\tau }_{b} $
      DRET10.69345.057.580.010.991.371.42%0.05%3.69%
      20.51366.2610.430.052.00%0.74%2.52%
      DGDL10.61343.897.16NA 10.972.6312.86%2.11%8.34%
      20.22361.2115.02NA 1NA 2NA 2NA 2
      30.31368.5510.98NA 118.97%2.31%6.22%
      1此处NA表示DGDL方法中不存在偏度系数。
      2此处NA表示该分量为错误分量,无法计算对应的参数误差。

      通过表2中的数据可以定量评估两种波形处理方法的分解效果:使用DRET方法分解的波形相关系数比DGDL方法提升了约2%,RMSE减小了48%。波形分量的幅值、峰值时刻和脉冲宽度的平均误差相较于DGDL方法分别降低了约14%、2%和4%。

      基于2000组已知参数的波形数据集,以进一步分析DRET波形处理算法的性能。具体处理结果的波形相关系数和RMSE的分布如图6所示。

      图  6  DRET和DGDL方法对已知参数波形数据集的处理结果

      Figure 6.  The processed results of the known-parameter waveform data set obtained by the DRET and DGDL methods

      图6中的所有处理结果按照分量间的峰值时刻间隔从大到小进行排列。图6显示,DRET方法所得到的波形相关系数和RMSE均优于DGDL方法。当波形分量的峰值时刻间隔较大时,DRET与DGDL方法都能对波形进行较好的分解,其波形相关系数均在较高水平,RMSE数值较小。但是,DGDL方法在间隔为8 ns左右出现了较大的分解误差,当分量峰值时刻间隔进一步减小时,DGDL方法的总体相关系数反而上升,这是由于两个波形分量的峰值时刻间隔很近时导致接收波形的形态更趋近于单分量波形,从而使得DGDL方法可以实现接收波形数据的准确拟合,却无法识别混叠的波形分量。DRET方法的总体性能稳定,相关系数和RMSE分布都较为集中,几乎不存在波动。

      波形特征参数的相对误差和分量提取正确率是体现波形处理算法性能的核心评价指标。统计不同峰值时刻间隔范围内平均的波形特征参数误差和分量提取正确率,其结果如图7所示。

      图  7  DRET和DGDL方法的波形特征参数提取误差和分量个数提取正确率

      Figure 7.  Extracted errors of the waveform parameters and the successful detection rates of the component numbers by the DRET and DGDL methods

      图7显示,DRET方法所得到的波形特征参数的相对误差明显小于DGDL方法。随着波形分量峰值时刻间隔的减小,由DGDL方法得到的波形特征参数相对误差呈现逐渐上升的趋势,而DRET方法的结果保持相对稳定。同时,DGDL方法提取得到分量正确率随着分量峰值时刻间隔的减小而大幅下降,而DRET方法得到的分量正确率始终保持在90%以上。统计两种方法所得到的评估指标参数的平均值见表3

      表 3  两种方法的评估指标平均值

      Table 3.  Average value of the evaluation indexes of the two methods

      Evaluation indexDRETDGDL
      $ {C}_{x} $ 0.987 0.966
      $ {\delta }_{x} $ 1.217 2.305
      $ {\tau }_{I} $ 2.18% 7.56%
      $ {\tau }_{u} $ 0.52% 3.69%
      $ {\tau }_{b} $ 2.33% 6.47%
      $ {r}_{s} $ 98.70% 64.25%

      表3中数据可知,DRET方法所得到的各项评估指标平均值都优于DGDL方法。其中,平均波形相关系数较DGDL方法提升了约2%,平均波形RMSE降低了约47%。DRET方法得到的所有波形特征参数相对误差平均值比DGDL方法降低了5%,波形分量提取正确率提高了约34%。

    • 在实验区域选取了1000组基于机载点云数据的仿真波形[25]作为分析对象,采用DRET和DGDL两种方法进行波形分解实验,得到的波形相关系数和RMSE的结果如图8所示。

      图  8  DRET方法和DGDL方法对机载仿真数据的分解结果

      Figure 8.  Decomposition results of airborne simulation data by the DRET and DGDL methods

      图8中,两种方法所得到的波形相关系数大部分集中在0.98以上,RMSE基本都小于4倍的噪声标准差。其中,DRET方法分解波形的相关系数和RMSE均优于DGDL方法,相关指标的稳定性更好。波形相关系数和RMSE的最小值、最大值、平均值和标准差的统计结果见表4

      表 4  DRET方法和DGDL方法对于机载仿真数据处理结果的相关系数($C_x $)和RMSE($\delta_x $)

      Table 4.  Correlation coefficient ($C_x $) and RMSE ($\delta_x $) of the DRET and the DGDL methods for the processing results of the airborne simulation data

      Evaluation indexDRETDGDL
      $ {C}_{x} $Max0.9990.991
      Min0.9370.872
      Average0.9950.986
      Standard deviation0.0120.034
      $ {\delta }_{x} $Max9.32614.273
      Min0.6071.988
      Average1.8693.516
      Standard deviation1.2912.936

      表4数据显示,DRET方法所得到的波形相关系数总体提升幅度不大,其均值仅比DGDL方法提升了约1%。但是两种方法所得到的RMSE指标差异明显,DRET方法的RMSE平均值相较于DGDL方法降低了约56%,说明DRET方法能够获得更好的波形分解效果,波形特征参数的准确性也更高。

    • 选取实验区域内GEDI轨迹上的1000个实际回波数据,采用DRET和DGDL方法进行波形处理,得到的波形相关系数和RMSE结果分布如图9所示。

      图  9  DRET和DGDL方法对GEDI数据的分解结果

      Figure 9.  Decomposition results of the GEDI data by the DRET and DGDL methods

      图9显示,GEDI数据波形分解的统计结果与机载仿真数据分解的统计结果分布类似,即DRET方法的波形相关系数和RMSE指标均优于DGDL方法,指标分布更加稳定。同理,统计得到波形相关系数和RMSE的最小值、最大值、平均值和标准差见表5

      表 5  DRET和DGDL方法对于GEDI数据处理结果的相关系数($C_x $)和RMSE ($\delta_x $)

      Table 5.  Correlation coefficient ($C_x $) and RMSE ($\delta_x $) of the DRET and the DGDL methods for the processing results of the GEDI data

      Evaluation indexDRETDGDL
      $ {C}_{x} $Max0.9990.982
      Min0.9390.854
      Average0.9930.977
      Standard deviation0.0370.062
      $ {\delta }_{x} $Max9.41215.504
      Min0.4493.075
      Average1.9534.248
      Standard deviation1.3633.712

      表5显示,DRET方法得到波形相关系数平均值仅比DGDL方法提升了约2%,而RMSE平均值比DGDL方法降低了54%。

    • 采用陡坡区域的植被机载点云数据、GEDI系统参数和波形仿真器,基于DRET方法所得到的波形特征参数结果提取植被高度用于定量评估该方法在目标参数反演方面的优势。植被目标的点云分布如图10所示。

      图  10  植被目标的点云分布

      Figure 10.  Point cloud distributions of the vegetated target

      图10中的被测目标由4棵树和地面组成,其中,地面目标的坡度设置为30°,最高树木的高度范围为5~30 m,点云中还存在其他低矮灌木。

      将机载点云数据中地面中心点与植被冠层顶部之间的高度作为真值,而由波形分解算法提取到的地面与植被波形分量的峰值时刻间隔作为植被高度解算值,得到不同树木高度条件下所解算的植被高度值与实际高度值差异的分布如图11所示。

      图  11  不同树高下植被高度解算的误差分布

      Figure 11.  Distributions of the canopy height errors under different tree heights

      图11显示,在不同树高条件下,DRET方法得到的植被高度误差均小于DGDL方法。其中,DGDL方法解算的植被高度误差最大值和平均值分别达到3.7 m和2.1 m,而DRET方法将最大误差和平均误差减小至1.2 m和0.6 m。

    • 文中提出了一种基于偏正态模型、激励的RL反卷积算法、逐层分解算法、梯度下降法和非线性最小二乘拟合相组合的星载激光测高仪接收波形的DRET分解方法。该方法特别适用于波形分量存在一定混叠以及形状偏离标准高斯形态而出现拖尾效应的多模式回波处理。

      基于已知波形参数的数据集、仿真波形数据集和GEDI实测波形数据的实验结果表明,波形相关系数和RMSE均优于传统的DGDL方法,其数值分布更为稳定;波形分量的特征参数的提取精度和分量提取正确率也明显高于DGDL方法。同时,基于DRET方法提取得到的植被高度也更为准确。

      DRET方法能够准确识别所有的波形分量并提取高精度的波形特征参数,这有利于该方法为线性体制星载激光测高仪的应用提供技术基础。

参考文献 (25)

目录

    /

    返回文章
    返回