CN110837088B - 一种星载激光测高仪数据去噪方法 - Google Patents

一种星载激光测高仪数据去噪方法 Download PDF

Info

Publication number
CN110837088B
CN110837088B CN201911054629.3A CN201911054629A CN110837088B CN 110837088 B CN110837088 B CN 110837088B CN 201911054629 A CN201911054629 A CN 201911054629A CN 110837088 B CN110837088 B CN 110837088B
Authority
CN
China
Prior art keywords
data
laser altimeter
intrinsic mode
component
satellite
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN201911054629.3A
Other languages
English (en)
Other versions
CN110837088A (zh
Inventor
谢欢
张志杰
童小华
刘世杰
许雄
陈鹏
金雁敏
冯永玖
柳思聪
王超
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Tongji University
Original Assignee
Tongji University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Tongji University filed Critical Tongji University
Priority to CN201911054629.3A priority Critical patent/CN110837088B/zh
Publication of CN110837088A publication Critical patent/CN110837088A/zh
Application granted granted Critical
Publication of CN110837088B publication Critical patent/CN110837088B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S17/00—Systems using the reflection or reradiation of electromagnetic waves other than radio waves, e.g. lidar systems
    • G01S17/02—Systems using the reflection of electromagnetic waves other than radio waves
    • G01S17/06—Systems determining position data of a target
    • G01S17/08—Systems determining position data of a target for measuring distance only
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2218/00—Aspects of pattern recognition specially adapted for signal processing
    • G06F2218/02—Preprocessing
    • G06F2218/04—Denoising

Landscapes

  • Physics & Mathematics (AREA)
  • Electromagnetism (AREA)
  • Engineering & Computer Science (AREA)
  • Computer Networks & Wireless Communication (AREA)
  • General Physics & Mathematics (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Remote Sensing (AREA)
  • Position Fixing By Use Of Radio Waves (AREA)

Abstract

本发明涉及一种星载激光测高仪数据去噪方法，其特征在于，该方法包括以下步骤：S1：获取星载激光测高仪数据；S2：基于经验模态分解，将星载激光测高仪数据分解为若干本征模态分量和一个单调余项；S3：计算每个本征模态分量的Hurst指数值；S4：基于每个本征模态分量的Hurst指数值，采用预建立的优化策略，获取去噪后的星载激光测高仪数据。与现有技术相比，本发明数据去噪方法能自适应地星载激光测高仪数据进行去噪，具有稳定性好、准确度高等优点。

Description

一种星载激光测高仪数据去噪方法
技术领域
本发明涉及数据处理领域，尤其是涉及一种星载激光测高仪数据去噪方法。
背景技术
星载激光测高仪数据处理中，激光全波形的处理和分解是提升测距精度和准确反演地物目标特征参数(如坡度、粗糙度等)的关键。由于受背景噪声、热噪声等的存在将对波形波形特征参数提取会造成干扰进而影响其测距精度，因此，激光测高全波形去噪方法的研究具有重要意义。目前，对全波形去噪主要方法有傅里叶低通滤波、小波滤波等，傅里叶低通滤波能有效消除高频噪声，但是能反映原始波形丰富细节程度的截至频率设定需要不断试验确定；不同的小波基函数、分解尺度、重构方法等去噪参数的选择对小波去噪去噪效果影响较大，且选定的小波基函数用于分析所有的波形数据时不具有自适应性。高斯滤波去噪效果依赖于核函数的宽度，过小的宽度去噪效果较弱，过大的宽度又会造成有效信号的损失。无法兼顾去噪效果和波形特征保留，对波形幅值的削弱和宽度展开会对波形分解的准确性产生影响；均值滤波不需要波形的先验知识，可以快速实现海量数据的处理，但也会造成回波数据的变形。
发明内容
本发明的目的就是为了克服上述现有技术存在对去噪参数依赖的缺陷而提供一种具有自适应性的星载激光测高仪数据去噪方法。
本发明的目的可以通过以下技术方案来实现：
一种星载激光测高仪数据去噪方法，该方法包括以下步骤：
S1：获取星载激光测高仪数据，该星载激光测高仪数据为时间序列；
S2：基于经验模态分解，将星载激光测高仪数据分解为若干本征模态分量和一个单调余项；
S3：计算每个本征模态分量的Hurst指数值；
S4：基于每个本征模态分量的Hurst指数值，采用预建立的优化策略，获取去噪后的星载激光测高仪数据。
进一步地，所述步骤S2具体包括以下步骤：
S201：将星载激光测高仪数据作为待分解信号；
S202：获取待分解信号的极大值点和极小值点；
S203：对所有极大值点进行拟合，获取上包络线，对所有极小值点进行拟合，获取下包络线；
S204：计算待分解信号的包络均值；
S205：从待分解信号中减去包络均值，从而对待分解信号进行更新；
S206：判断更新后的待分解信号是否满足本征模态分量条件，若不满足，则重复步骤S202至S205，若满足，则获得一个本征模态分量；
S207：从星载激光测高仪数据中减去步骤S206获取的本征模态分量，从而对星载激光测高仪数据进行更新，重复步骤S201至S206，对更新后的星载激光测高仪数据迭代进行经验模态分解，获取若干本征模态分量和一个单调余项。
进一步地，所述步骤S203中，采用spline插值方法获取上包络线和下包络线。
进一步地，所述步骤S206中，判断更新后的待分解信号是否满足本征模态分量条件具体为，若更新后的待分解信号相邻两次迭代之间的标准差为0.25，则满足本征模态分量条件，否则不满足本征模态分量条件，所述标准差的计算表达式为：
Figure BDA0002256222670000021
式中，SD为标准差，T为一本征模态分量中时间节点的个数，IMF1(k)(t)为该本征模态分量第k次迭代时，第t个时间节点的值，IMF1(k-1)(t)为该本征模态分量第k-1次迭代时，第t个时间节点的值。
进一步地，所述步骤S3具体为，采用拟合阶数为1的消除趋势波动分析法，计算每个本征模态分量的Hurst指数值。
进一步地，所述拟合阶数为1的消除趋势波动分析法具体包括以下步骤：
S301：计算某一本征模态分量的累计离差，滤除本征模态分量的平均值，对本征模态分量进行更新；
S302：将更新后的本征模态分量划分为不相交的等长的子区间；
S303：计算每个子区间的局部趋势，通过将每个子区间滤除局部趋势，对每个子区间进行更新；
S304：计算更新后的每个子区间的方差；
S305：基于方差，对所有子区间求均值并开方，获取整个星载激光测高仪数据的波动函数；
S306：基于波动函数，获取本征模态分量的Hurst指数值；
S307，选取另一本征模态分量，重复步骤S301至S306，直至遍历所有本征模态分量。
进一步地，所述步骤S4中，优化策略具体为，若本征模态分量的Hurst指数值小于或等于预建立的Hurst指数阈值，则去除该本征模态分量，否则保留该本征模态分量。
进一步地，所述Hurst指数阈值为0.5。如果Hurst指数小于0.5，则表明本征模态分量的时间序列的总体趋势与过去相反；如果Hurst指数等于0.5，则表明本征模态分量的时间序列是随机的，不会影响未来；如果Hurst指数大于0.5，则表明本征模态分量的时间序列具有持续性，Hurst指数越大，相关性、持续性就越强，因此将Hurst指数阈值设为0.5，即确保本征模态分量具有持续性，不存在反趋势，实现对数据的去噪。
与现有技术相比，本发明具有以下优点：
(1)本发明星载激光测高仪数据去噪方法，采用基于Huang等提出的经验模态分解方法进行数据分解，经验模态分解方法可以将复杂信号由精细尺度到粗大尺度分解为若干本征模态分量和一个单调余项，通常全波形中噪声的频率远大于回波高斯脉冲的频率，因此可以通过舍弃一定高频的本征模态分量实现去噪的目的，本发明通过对每个本征模态分量计算Hurst指数实现高频噪声的判断，该方法考虑了每个本征模态分量本身的特点实现了本征模态分量选择的自适应性，避免了对滤波参数的依赖。
(2)仿真试验结果表明，本发明基于经验模态分解和Hurst指数的星载激光测高仪数据去噪方法，可以全自动筛选本征模态分量实现去噪目的，同时对不同信噪比波形具备较好的稳定性。对ICESat/GLAS真实数据拟合效果更优，且峰值个数判断更为精确。
(3)本发明星载激光测高仪数据去噪方法，考虑到采用消除趋势波动分析法进行Hurst指数的计算时，拟合阶数K越大会对较短时间序列Hurst指数估算引入较大的误差，因此将拟合阶数设为1，确保本征模态分量的Hurst指数计算结果的准确性。
附图说明
图1为本发明星载激光测高仪数据去噪方法的流程示意图；
图2为本实施例中不同筛选IMFs方法去噪效果对比图，其中，(a)为不同筛选IMFs方法去噪前的幅值强度图，(b)为不同筛选IMFs方法去噪后的幅值强度图；
图3为本实施例中信噪比增加后波形分解一致率变化图；
图4为本实施例中信噪比增加后波形分解一致率分布图；
图5为本实施例中不同方法的第一拟合参数误差分布图，其中每种方法对应的四个长条形框从左往右依次为A、T、σ和FitDev的数据值；
图6为本实施例中不同方法的第二拟合参数误差分布图，其中每种方法对应的四个长条形框从左往右依次为A、T、σ和FitDev的数据值；
图7为本实施例中一条典型GLAH01波形所有的IMFs及余项图；
图8为本实施例中分解平均一致率与GLAH05产品参数对比图；
图9为本实施例中EMD-Hurst fitting与GLAH05高斯成分个数对比图；
图10为本实施例中EMD-Hurst fitting与GLAH05第一拟合结果对比图，其中(a)为GLAH05拟合结果图，(b)为EMD-Hurst fitting拟合结果图；
图11为本实施例中EMD-Hurst fitting与GLAH05第二拟合结果对比图，其中(a)为GLAH05拟合结果图，(b)为EMD-Hurst fitting拟合结果图。
具体实施方式
下面结合附图和具体实施例对本发明进行详细说明。本实施例以本发明技术方案为前提进行实施，给出了详细的实施方式和具体的操作过程，但本发明的保护范围不限于下述的实施例。
实施例1
目前针对星载激光测高全波形数据去噪方法主要有固定宽度的高斯滤波器、傅里叶低通滤波、小波滤波等，以上滤波器的效果均依赖于滤波参数的选择且不具有自适应性。
针对这个问题，本实施例一种星载激光测高仪数据去噪方法基于具有自适应经验模态(Empirical Mode Decomposition，EMD)分解方法，对分解后的每个本征模态分量(IMF)进行Hurst(H)指数分析，将H≤0.5的IMF归为高频噪声成分给与直接去除，该方法相比传统基于EMD分解去噪方法考虑了每个IMF本身的特点实现了IMF选择的自适应性。为了验证该方法的有效性，将其去噪效果与其他IMF筛选方法如阈值函数(EDM-Threshold)、相关性法(EMD-Correlation)、结合小波滤波器法(EMD-Wavelets)、直接去除前N个IMFs等6种方法进行对比，并利用大批量仿真波形和ICESat/GLAS波形数据进行了验证与分析。仿真实验结果表明，EMD-Hurst指数可以全自动筛选IMFs实现去噪目的，同时对不同信噪比波形具备较好的稳定性。对ICESat/GLAS真实数据拟合效果更优，且峰值个数判断更为精确。
本实施例一种星载激光测高仪数据去噪方法，包括以下步骤：
S1：获取星载激光测高仪数据；
S2：基于经验模态分解，将星载激光测高仪数据分解为若干本征模态分量和一个单调余项；
S3：计算每个本征模态分量的Hurst指数值；
S4：基于每个本征模态分量的Hurst指数值，采用预建立的优化策略，获取去噪后的星载激光测高仪数据。
下面对本实施例一种星载激光测高仪数据去噪方法进行具体介绍：
1、经验模态分解
EMD分解最早由Huang提出对非线性、非平稳信号进行分析，该方法认为任何复杂的时间序列都是若干阶不相同的、简单的、非正弦函数的IMF组成，基于此可以从复杂时间序列直接分离出从高频到低频基本时间序列的IMFs。每个IMF需要满足两个条件：(1)整个IMF中零点数与极点数相等或至多相差1；(2)IMF上任意一点由局部极大值点确定的上包络线和由局部极小值点确定的下包络线均值都为零。获取每个IMF成分的筛选过程需要迭代满足以上两个条件，然后反复从原始全波形数据中减去已分离的IMFs直至最后一个余项呈现单调趋势，即完成了全波形EMD的整个分解过程。在本实施例研究中，采用spline插值方法对上下包络线进行插值，分解过程终止条件相邻两次迭代之间标准差取0.25。
经验模态分解的过程具体包括以下步骤：
S201：将星载激光测高仪数据wf(t)作为待分解信号s(t)＝wf(t)；
S202：获取待分解信号s(t)的极大值点和极小值点；
S203：，采用spline插值方法对所有极大值点进行插值拟合，获取上包络线h(t)，采用spline插值方法对所有极小值点进行插值拟合，获取下包络线l(t)；
S204：计算待分解信号的包络均值m(t)，m(t)＝(h(t)+l(t))/2；
S205：从待分解信号中减去包络均值，从而对待分解信号进行更新；
对待分解信号更新过程的计算公式为c11(t)＝s(t)-m(t)
S206：判断更新后的待分解信号c11(t)是否满足本征模态分量条件，若不满足，则重复步骤S202至S205，求c12(t)、c13(t)、…、c1k(t)，若满足，则获得第一个本征模态分量IMF1；
判断更新后的待分解信号是否满足本征模态分量条件具体为，若更新后的待分解信号相邻两次迭代之间的标准差为0.25，则满足本征模态分量条件，否则不满足本征模态分量条件，标准差的计算表达式为：
Figure BDA0002256222670000061
式中，SD为标准差，T为一本征模态分量中时间节点的个数，IMF1(k)(t)为该本征模态分量第k次迭代时，第t个时间节点的值，IMF1(k-1)(t)为该本征模态分量第k-1次迭代时，第t个时间节点的值。
S207：从星载激光测高仪数据中减去步骤S206获取的本征模态分量，从而对星载激光测高仪数据进行更新s(t)＝s(t)-IMFi，IMFi为第i个本征模态分量，重复步骤S201至S206，对更新后的星载激光测高仪数据迭代进行经验模态分解，获取若干本征模态分量和一个单调余项r(t)。
2、Hurst指数
Hurst指数分析法通常用来分析时间序列的分形特征和长期记忆过程，最初是由英国水文专家H.E.Hurst在研究尼罗河水库水流量和储存能力的关系时提出，而后它被用在各种时间序列的分析之中。目前，国内外学者已经提出了多种估计Hurst指数的方法。其中Hurst提出了重极标差(R/S)分析法，以判断时间序列数据遵从随机游走是有偏随机游走过程的指标。Peng等在研究DNA组织时扩展了普通的波动分析方法，得到了(DetrendedFluctuation Analysis，DFA)方法，它在消除时间序列局部趋势及发现局部相关方面比R/S及改进R/S方法更优。Kantelhardt等对DFA拟合阶数DFA-K分析得到拟合阶数K越大会对较短时间序列Hurst指数估算引入较大的误差。故本实施例采用DFA-1，即线性回归分析来构建Hurst指数。设获取的星载激光测高仪数据为长度为N的时间序列{x_t，t＝1，2，3…N}，对时间序列进行经验模态分解，得到若干本征模态分量，从而计算每个本征模态分量的Hurst指数。
如图1所示，Hurst指数的计算包括以下步骤：
S301：计算某一本征模态分量的累计离差，滤除本征模态分量的平均值，对本征模态分量进行更新，更新过程的计算表达式为
Figure BDA0002256222670000071
式中，Y(i)为更新后的本征模态分量，
Figure BDA0002256222670000072
为本征模态分量总值，
Figure BDA0002256222670000073
为本征模态分量平均值。
S302：将更新后的本征模态分量Y(i)划分为不相交的长度为s的Ns＝int(N/s)个等长的子区间；
S303：计算每个子区间的局部趋势，通过将每个子区间滤除局部趋势，对每个子区间进行更新；
局部趋势的计算具体为，用最小二乘法对子区间所包含的数据进行一阶线性拟合，获取局部趋势Ys；
S304：计算更新后的每个子区间的方差；
方差计算的表达式为
Figure BDA0002256222670000074
式中，
Figure BDA0002256222670000075
为长度为s的第v段子区间的方差，v＝1,…,Ns，Ns为子区间的总数；
S305：基于方差，将所有子区间求均值并开方，获取整个星载激光测高仪数据的波动函数；
S306：基于波动函数，获取本征模态分量的Hurst指数值；
S307，选取另一本征模态分量，重复步骤S301至S306，直至遍历所有本征模态分量。
3、试验方案
为了验证提出方法的有效性，本实施例采用仿真和真实的波形数据进行实验。其中仿真波形中，先定性观察EMD-Hurst(本实施例星载激光测高仪数据去噪方法)和其他6种基于EMD分解方法在不同信噪比条件下去噪效果，然后对批量波形通过分解后的评价指标统计量进行定量地评价EMD-Hurst方法在不同信噪比条件下稳健性。同样，在海冰厚度，冰、云和陆地高程卫星/地球科学激光测高系统(ice，cloud and land elevation satellite/geoscience laser altimeter system，ICESat/GLAS)获取的真实波形中，首先结合EMD分解后的IMF及其相应的Hurst指数，观察分析了去除IMF是否高频噪声部分。然后批量地对比不同的方法的分解情况，最后将基于EMD-Hurst指数去噪波形分解结果与GLAH05产品结果进行对比分析，统计一致率并分析不一致率原因。
3.1、去噪评价指标
传统评价指标主要有信噪比和去噪前后波形的均方根误差，有这两个指标的定义可知并不是信噪比越大越好，均方根误差越小越好。故本实施例基于发射脉冲符合高斯分布，假设单个地物回波也符合高斯分布，把整个回波看作是若干个高斯函数的叠加，对滤波后回波按照高斯模型(Levenberg-Marquardt，LM)拟合分解后，提出一些客观定量评价指标：①对仿真波形参数，设置分解一致率，其满足两个条件：分解波形个数一致且峰值位置误差不超过1ns；鉴于真值已知，另外设置峰值幅值A、峰值位置T、峰值宽度σ均值误差对比。②对ICESat/GLAS波形，通过正常波形个数比总波形个数设置总分解一致率、分解后波峰个数相对于总回波个数设置波峰个数一致率、峰值个数一致率比总分解一致率设置相对分解一致率。
3.2、仿真试验设计
本实施例的仿真数据是在matlab平台上，在仿真混合高斯模型加上加上随机高斯白噪声。信噪比从10db以步长5间隔递增至40db。回波采样个数为544，样本间隔1ns。最大峰值个数3个，幅值强度(Amp)大小1-10，位置(Loc)在100-500ns范围内分布，宽度(Sigma)最小值为6ns。具体信噪比(SNR)、幅值水平、峰值位置、峰值宽度设置参数如表1所示。为了避免少量回波结果的偶然性，通过变换峰值参数和信噪比生成2000条仿真波形并对实验结果进行了统计分析。
表1仿真波形参数设置
Figure BDA0002256222670000091
3.2.1、仿真波形去噪比较
随机产生一条波形，通过对其建立Hurst指数、相关性分析、直接去除第一个IMF、直接去除前2个IMFs、去除前两个IMFs加上软硬阈值、对分解后的所有的IMFs进行小波滤波，其中结合星载激光测高全波形的特点，本实施例选取dn4为小波基，对EMD分解之后的所有IMFs进行3层分解，软阈值法进行信号重构，即EMD-Wavelet(小波滤波器法)。各种筛选IMFs方法去噪效果对比如图2所示，由图2可知(1)去噪效果最明显是EMD-Hurst指数、EMD-Wavelet；(2)直接去除第一个IMF去噪效果不足，相比直接去除第一个IMF，直接去除前两个IMF或再加上软硬阈值方法均取得了较好的去噪效果；(3)去噪效果最差是基于相关性分析，原因可能是硬阈值设置不具备很好的适应性。
3.2.2、仿真试验结果与讨论
为了验证各种筛选IMFs方法去噪效果，随着信噪比增加，几种筛选IMFs方法波形分解一致率变化情况如图3所示。由图3可知，(1)相比其他方法，基于EMD-Hurst指数分析方法，在不同信噪比条件下均取得了较高的分解一致率。(2)基于EMD-Correlation方法随着信噪比增加不具备稳定性，不同信噪比条件下，分解一致率相差较大；(3)直接去除第一个IMF，因去噪效果不足，整个分解一致率较低。(4)直接去掉前2个IMFs或者加上软硬阈值进行修正后，分解一致率优于直接去除第一个IMF。(5)除了EMD-Wavelet随着信噪增加分解一致率一直上升外，其余各种方法均出现了分解一致率先上升后下降的趋势。
由各种筛选方法实现的原理可分析得，这七种方法可以归为三类：(1)能对EMD分解后IMFs实现全自动判别的EMD-Hurst、EMD-Correlation(相关性法)；(2)对IMFs半自动判别的EMD-1stIMF(EMD-1IMF，去除第一个IMF法)、EMD-First2IMFs(EMD-2IMFs，去除第二个IMF法)、EMD-soft(EMD-softThresholding，软阈值法)、EMD-hard(EMD-hardThresholding，硬阈值法)；(3)对IMFs不判别的EMD-Wavelet。其中，基于第一类筛选方法，EMD-Correlation在不同信噪比波形最大的相关性不一致，取η＝5设置硬阈值λ不能适应波形复杂度和信噪比变化；EMD-Hurst指数分析，则能较好且稳定地适应波形信噪比变化。第二类筛选方法，依赖于前N个高频IMFs设定，N设置过小则去噪效果不足，N过大则会扔掉有效信号；加上软硬阈值之后会使得直接去掉的有效信号得到一定程度补偿；第三类，无论信噪比高低都进行去噪，故其分解一致率随着信噪比增加在提高。相较于第三类，第一类方法则随着信噪比提高，筛选的IMFs的Hurst指数或硬阈值λ均大于设定阈值，没有去除IMFs即没有实现对原波形去噪平滑，故分解一致率会随着信噪比提高而降低。第二类随着信噪比增加，由于强制去除前N个IMFs一定会引起有效信号的损失，故分解一致率也呈现随着信噪比增加而下降的趋势。
如图4至6所示，为了综合对比各种筛选方法的效果，本实施例计算了上述评价指标分解一致率CR大小、峰值三个参数A、T、σ绝对误差、拟合均方根误差FitDev的均值，由图4至6可知，EMD-Hurst指数分析，平均分解一致率最高；除了EMD-Correlation、EMD-1stIMF分解一致率较低，波峰参数平均误差较大外，其余方法误差相当。
3.3、ICESat/GLAS波形
为了使结果更有说服力，同时验证仿真结果。本实施例随机选择了1000条GLAH01回波波形和相应的GLAH05产品分解参数进行验证分析。其中，在波形分解过程种，基于ICESat/GLAS系统实际情况，对分解后的波峰提出如下约束条件：(1)峰值最小值大于噪声均值加3倍噪声中误差；(2)峰值位置分布在有效信号区间内；(3)幅值半宽范围(2.5<σ<300)；(4)相邻峰值位置相差大于15ns。
3.3.1、EMD-Hurst去噪效果
一条典型ICESat/GLAS波形分解后的6个IMFs如图7所示，其对应的Hurst指数分别是0.2838，0.3145，0.6988，0.8826，0.9645，1.0140。由图7可知高频部分主要集中在前两个IMFs，这与Hurst指数评价结果相一致。
3.3.2、波形分解结果与讨论
如图8所示，基于上述随机选择的1000条GLAH01回波波形，将不同基于EMD去噪方法后拟合结果和相应的GLAH05产品分解参数进行验证分析。由图8可知，各种去除IMFs方法总分解不一致率均没有到达1，是源于σ小于分解设定阈值3ns，被LM优化掉。
3.3.3、EMD-Hurst与GLAH05对比结果与讨论
如图9至对比基于EMD-Hurst和GLAH05产品参数对原始波形进行拟合对比发现引起原因GLAH05参数拟合效果较差且分解参数不合理。如图10所示，EMD-Hurst波峰2个，GLAH05波峰4个；如图11所示，EMD-Hurst波峰3个，GLAH05波峰5个；可知GLAH05拟合效果偏差较大。对应参数表2和表3表明GLAH05过多估计峰值个数且峰值参数存在不合理数值(加粗参数)，如相邻峰值位置小于系统测距分辨率1ns，系统测距分辨率即可分辨率1ns/15cm高差。这种情况下，如果基于最后一个或2个均值峰值反演确定地面位置势必会造成较大的误差。
表2.EMD-Hurst与GLAH05产品峰值参数对比(对应图10)
Figure BDA0002256222670000111
表3.EMD-Hurst与GLAH05产品峰值参数对比(对应图11)
Figure BDA0002256222670000112
4、结论
本实施例基于EMD分解去噪原理，针对分解后高频IMFs筛选问题借鉴其他领域对非稳态信号的处理方法，提出构建Hurst指数方法来识别该IMF是否是需要去除的高频噪声。通过大量仿真波形实验可分析：(1)对比其他半自动设置直接去除前N个IMF或者加阈值函数进行补偿方法，该方法避免了由于N设置过小去噪不足，N过大造成有效信号损失问题。(2)随着信噪增加，波峰分解一致率均值最大、综合去噪效果最稳定。此外，与GLAH05产品参数对比分析，基于EMD-Hurst指数分析法波形拟合精度更高且峰值分解个数更准确。根据本实施例介绍的方法，可以应用到不同信噪比的星载激光测高波形去噪中，为星载激光测高波形去噪处理提出了新的处理方法。
以上详细描述了本发明的较佳具体实施例。应当理解，本领域的普通技术人员无需创造性劳动就可以根据本发明的构思作出诸多修改和变化。因此，凡本技术领域中技术人员依本发明的构思在现有技术的基础上通过逻辑分析、推理或者有限的实验可以得到的技术方案，皆应在由权利要求书所确定的保护范围内。

Claims (7)

1.一种星载激光测高仪数据去噪方法，其特征在于，该方法包括以下步骤：
S1：获取星载激光测高仪数据；
S2：基于经验模态分解，将星载激光测高仪数据分解为若干本征模态分量和一个单调余项；
S3：计算每个本征模态分量的Hurst指数值；
S4：基于每个本征模态分量的Hurst指数值，采用预建立的优化策略，获取去噪后的星载激光测高仪数据；
所述步骤S2具体包括以下步骤：
S201：将星载激光测高仪数据作为待分解信号；
S202：获取待分解信号的极大值点和极小值点；
S203：对所有极大值点进行拟合，获取上包络线，对所有极小值点进行拟合，获取下包络线；
S204：计算待分解信号的包络均值；
S205：从待分解信号中减去包络均值，从而对待分解信号进行更新；
S206：判断更新后的待分解信号是否满足本征模态分量条件，若不满足，则重复步骤S202至S205，若满足，则获得一个本征模态分量；
S207：从星载激光测高仪数据中减去步骤S206获取的本征模态分量，从而对星载激光测高仪数据进行更新，重复步骤S201至S206，对更新后的星载激光测高仪数据迭代进行经验模态分解，获取若干本征模态分量和一个单调余项；
所述步骤S206中，判断更新后的待分解信号是否满足本征模态分量条件具体为，若更新后的待分解信号相邻两次迭代之间的标准差为0.25，则满足本征模态分量条件，否则不满足本征模态分量条件，所述标准差的计算表达式为：
Figure FDA0003953648780000011
式中，SD为标准差，T为一本征模态分量中时间节点的个数，IMF1(k)(t)为该本征模态分量第k次迭代时，第t个时间节点的值，IMF1(k-1)(t)为该本征模态分量第k-1次迭代时，第t个时间节点的值。
2.根据权利要求1所述的一种星载激光测高仪数据去噪方法，其特征在于，所述步骤S203中，采用spline插值方法获取上包络线和下包络线。
3.根据权利要求1所述的一种星载激光测高仪数据去噪方法，其特征在于，所述步骤S3具体为，采用拟合阶数为1的消除趋势波动分析法，计算每个本征模态分量的Hurst指数值。
4.根据权利要求3所述的一种星载激光测高仪数据去噪方法，其特征在于，所述拟合阶数为1的消除趋势波动分析法具体包括以下步骤：
S301：计算某一本征模态分量的累计离差，滤除本征模态分量的平均值，对本征模态分量进行更新；
S302：将更新后的本征模态分量划分为不相交的等长的子区间；
S303：计算每个子区间的局部趋势，通过将每个子区间滤除局部趋势，对每个子区间进行更新；
S304：计算更新后的每个子区间的方差；
S305：基于方差，对所有子区间求均值并开方，获取整个星载激光测高仪数据的波动函数；
S306：基于波动函数，获取本征模态分量的Hurst指数值；
S307，选取另一本征模态分量，重复步骤S301至S306，直至遍历所有本征模态分量。
5.根据权利要求4所述的一种星载激光测高仪数据去噪方法，其特征在于，所述局部趋势的计算过程具体为，用最小二乘法对子区间所包含的数据进行一阶线性拟合，获取局部趋势。
6.根据权利要求1所述的一种星载激光测高仪数据去噪方法，其特征在于，所述步骤S4中，优化策略具体为，若本征模态分量的Hurst指数值小于或等于预建立的Hurst指数阈值，则去除该本征模态分量，否则保留该本征模态分量。
7.根据权利要求6所述的一种星载激光测高仪数据去噪方法，其特征在于，所述Hurst指数阈值为0.5。
CN201911054629.3A 2019-10-31 2019-10-31 一种星载激光测高仪数据去噪方法 Active CN110837088B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201911054629.3A CN110837088B (zh) 2019-10-31 2019-10-31 一种星载激光测高仪数据去噪方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201911054629.3A CN110837088B (zh) 2019-10-31 2019-10-31 一种星载激光测高仪数据去噪方法

Publications (2)

Publication Number Publication Date
CN110837088A CN110837088A (zh) 2020-02-25
CN110837088B true CN110837088B (zh) 2023-02-10

Family

ID=69575908

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201911054629.3A Active CN110837088B (zh) 2019-10-31 2019-10-31 一种星载激光测高仪数据去噪方法

Country Status (1)

Country Link
CN (1) CN110837088B (zh)

Families Citing this family (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111665517B (zh) * 2020-05-29 2022-11-18 同济大学 一种基于密度统计的单光子激光测高数据去噪方法及装置
CN111948663A (zh) * 2020-08-18 2020-11-17 中国科学院上海天文台 一种星载全波形信号的自适应经验模态分解去噪方法

Family Cites Families (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101833087B (zh) * 2010-05-13 2012-10-24 中国科学院对地观测与数字地球科学中心 基于小波分析的星载激光雷达波形数据分解的方法
CN103792526B (zh) * 2014-02-24 2016-03-30 武汉大学 一种基于脉冲回波形态的激光测高仪动态阈值选取方法
CN104635223A (zh) * 2015-01-08 2015-05-20 南京大学 基于经验模式分解与分数阶傅里叶变换的激光回波去噪方法
WO2017205382A1 (en) * 2016-05-23 2017-11-30 The University Of New Hampshire Techniques for empirical mode decomposition (emd)-based signal de-noising using statistical properties of intrinsic mode functions (imfs)

Also Published As

Publication number Publication date
CN110837088A (zh) 2020-02-25

Similar Documents

Publication Publication Date Title
CN101727662B (zh) Sar图像非局部均值去斑方法
Montanari et al. Estimating long-range dependence in the presence of periodicity: an empirical study
CN110837088B (zh) 一种星载激光测高仪数据去噪方法
CN110673109B (zh) 一种星载大光斑激光雷达全波形数据分解方法
CN111008585B (zh) 基于自适应分层高分辨sar图像的舰船目标检测方法
CN109871733A (zh) 一种自适应海杂波信号去噪方法
CN114580260B (zh) 一种基于机器学习和概率理论的滑坡区间预测方法
CN111679328A (zh) 一种基于变分模态分解的瞬变电磁探测信号的降噪方法
CN110988834A (zh) 一种基于自适应阈值双门限的脉冲到达时间测量方法
CN113917490A (zh) 激光测风雷达信号去噪方法及装置
CN108489529A (zh) 一种基于高阶统计量的微弱信号检测方法
Guo et al. Groundwater depth forecasting using configurational entropy spectral analyses with the optimal input
CN106526545A (zh) 一种鲁棒cfar检测器的检测方法
CN111830481B (zh) 雷达回波单分量幅度分布模型参数估计方法及装置
CN106772573B (zh) 基于最大相关熵的地震子波信号提取方法
CN116032709B (zh) 无先验知识fsk信号盲解调和调制特征解析方法及装置
Stadnytska et al. Analyzing fractal dynamics employing R
Kazakov et al. Algorithm of Robust Frequency Estimation in a Channel with White Gaussian Noise and Pulse Interferences
CN111077532A (zh) 一种基于反卷积和高斯分解的地物空间信息获取方法
CN115905919A (zh) 基于afgan的雷达信号样本数据处理方法、目标识别方法及系统
CN113052053B (zh) 一种坐标时间序列周期项的提取方法
CN112539142B (zh) 一种针对海上风电结构状态监测数据噪声的分析识别方法
CN106446548A (zh) 一种基于小波分析的动力学结构突变检测方法
CN103854658A (zh) 基于短时分形维数增强法的钢板腐蚀声发射信号去噪方法
CN114840802B (zh) 一种判别水文气候过程自然演变类型的方法

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant