CN107121705B - 一种基于自动反相校正和峰度值比较的探地雷达回波信号去噪算法 - Google Patents
一种基于自动反相校正和峰度值比较的探地雷达回波信号去噪算法 Download PDFInfo
- Publication number
- CN107121705B CN107121705B CN201710296899.XA CN201710296899A CN107121705B CN 107121705 B CN107121705 B CN 107121705B CN 201710296899 A CN201710296899 A CN 201710296899A CN 107121705 B CN107121705 B CN 107121705B
- Authority
- CN
- China
- Prior art keywords
- signal
- kurtosis value
- kurtosis
- value
- signals
- 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
Links
- 238000004422 calculation algorithm Methods 0.000 title claims abstract description 42
- 230000000149 penetrating effect Effects 0.000 title claims abstract description 15
- 238000012937 correction Methods 0.000 title claims abstract description 13
- 238000000034 method Methods 0.000 claims abstract description 9
- 238000012880 independent component analysis Methods 0.000 claims description 15
- 239000011159 matrix material Substances 0.000 claims description 11
- 238000004364 calculation method Methods 0.000 claims description 5
- 230000000694 effects Effects 0.000 abstract description 6
- 238000004458 analytical method Methods 0.000 abstract description 2
- 238000000354 decomposition reaction Methods 0.000 description 11
- 238000012545 processing Methods 0.000 description 6
- 238000011160 research Methods 0.000 description 4
- 238000001514 detection method Methods 0.000 description 3
- 238000002592 echocardiography Methods 0.000 description 3
- 239000002184 metal Substances 0.000 description 3
- 238000003672 processing method Methods 0.000 description 3
- 238000004088 simulation Methods 0.000 description 3
- 230000002159 abnormal effect Effects 0.000 description 2
- 239000000284 extract Substances 0.000 description 2
- 238000000926 separation method Methods 0.000 description 2
- 101100136092 Drosophila melanogaster peng gene Proteins 0.000 description 1
- 238000009933 burial Methods 0.000 description 1
- 238000010276 construction Methods 0.000 description 1
- 238000013461 design Methods 0.000 description 1
- 238000010586 diagram Methods 0.000 description 1
- 239000006185 dispersion Substances 0.000 description 1
- 238000000605 extraction Methods 0.000 description 1
- 238000005070 sampling Methods 0.000 description 1
- 238000012216 screening Methods 0.000 description 1
- 230000003595 spectral effect Effects 0.000 description 1
- 230000009466 transformation Effects 0.000 description 1
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V3/00—Electric or magnetic prospecting or detecting; Measuring magnetic field characteristics of the earth, e.g. declination, deviation
- G01V3/12—Electric or magnetic prospecting or detecting; Measuring magnetic field characteristics of the earth, e.g. declination, deviation operating with electromagnetic waves
-
- 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
- G01S7/00—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00
- G01S7/02—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S13/00
- G01S7/28—Details of pulse systems
- G01S7/285—Receivers
- G01S7/292—Extracting wanted echo-signals
Landscapes
- Engineering & Computer Science (AREA)
- Remote Sensing (AREA)
- Physics & Mathematics (AREA)
- Life Sciences & Earth Sciences (AREA)
- General Physics & Mathematics (AREA)
- Electromagnetism (AREA)
- Radar, Positioning & Navigation (AREA)
- Computer Networks & Wireless Communication (AREA)
- Environmental & Geological Engineering (AREA)
- Geology (AREA)
- General Life Sciences & Earth Sciences (AREA)
- Geophysics (AREA)
- Geophysics And Detection Of Objects (AREA)
Abstract
本发明公开了一种基于自动反相校正和峰度值比较的探地雷达回波信号去噪算法,包括以下步骤:将探地雷达原始回波信号与相同长度的零均值白噪声信号进行随机拟合,得到两道信号;运用独立成分分析算法对这两道信号进行分离处理,输出峰度值较大的信号记为x′(m),峰度值较小的信号记为n′(m);对x′(m)进行自动反相校正,再运用完全总体经验模态算法对反相校正后的信号进行分解得到P个信号分量,计算各信号分量的峰度值;再计算信号n′(m)的峰度值作为阈值;最后将峰度值大于阈值的所有信号分量累加,作为原始回波信号x(m)去噪后的信号。本方法提高了计算效率和去噪效果。
Description
技术领域
本发明属于探地雷达探测与应用技术领域,具体涉及探地雷达回波信号去噪处理。
背景技术
探地雷达(Ground Penetrating Rada,GPR)通过发射天线向地下发射宽带电磁波,接收天线接收散射回波,通过对散射回波进行处理,对地下未知区域进行无损探测和参数反演。电磁波在地下介质中传播时,遇到存在电性差异的分界面时发生散射,根据接收到的电磁散射回波来反演地下未知区域的异常***置、形态、埋藏深度等参数。GPR作为探测地下目标的一种重要工具,在人道主义救援中探测地下生命迹象,在维和行动中定位清除地雷、未***武器,在安全保卫中探测判断地下异常物,在施工中定位地下管道、线缆等方面正得到越来越多的研究与使用。
GPR回波信号由直达波、界面反射波、目标散射波、随机噪声等成分构成,在频率域和时间域互有重叠,难以区分。GPR发送的电磁波在地下传播的过程中,受到地下介质结构复杂多变、仪器参数或者噪声的多种因素的影响,会发生衰减、频散和其他干扰,这在很大程度上限制了GPR的探测分辨率和数据解释效果。因此,为了获得质量良好的回波信号,需要抑制噪声和杂波,提取所需的目标回波信号,最大限度地还原地下的异常体的散射回波。GPR信号处理水平对定位和识别目标起着决定性作用,因此开展GPR信号处理方法的研究工作有着重要的研究价值与现实意义。
国内外开展了大量的GPR信号去噪处理方面的研究工作。文献1“Jing Li,CaiLiu,Zhaofa Zeng,Lingna Chen,GPR Signal Denoising and Target Extraction Withthe CEEMD Method,IEEE Geoscience and Remote Sensing Letters,2015,12(8):1615-1619.”提出了一种基于完全总体经验模态分解(Complete Ensemble Empirical ModeDecomposition,CEEMD)的GPR信号处理方法,通过对GPR回波信号进行CEEMD分解,人工方式进行噪声辨识和去除,进而实现信号去噪。运用Hilbert-Huang变换,证明了该算法比经验模态分解和总体经验模态分解的谱分辨率更高;文献2“谢朋成,探地雷达探测浅地层目标的算法研究,东北林业大学,2016”提出了一种基于独立成分分析(Independent ComponentAnalysis,ICA)的GPR探测浅地层目标的算法研究。文献3“陈玲娜,基于CEEMD和PCA的探地雷达数据处理研究与应用,吉林大学,2016”对复杂GPR信号进行分频处理,通过选择固定频段的本征模态函数(Intrinsic Mode Function,IMF)进行重构从而去除高频杂波和低频隐失波的影响,再通过从原始信号中去掉高阶主成分信息,提取目标信号。上述的处理方法中,基于ICA的去噪算法存在信号相位不定性的问题,基于CEEMD的去噪算法则需要通过人工判读的方式对各IMF分量进行判别,再进行信号重构,效率低。
因此,有必要设计一种能够自动判定相位并自动进行IMF分量提取的GPR信号去噪算法。
发明内容
本发明所要解决的技术问题是,针对ICA算法的信号相位不定性和CEEMD算法中问题,提供一种基于自动反相校正和峰度值比较的探地雷达回波信号去噪算法,提高了计算效率和去噪效果。
本发明的技术方案如下:
一种基于自动反相校正和峰度值比较的探地雷达回波信号去噪算法,包括以下步骤:
步骤1:首先将探地雷达单道含噪原始回波信号x(m),m=1,2,...,M与相同长度的零均值白噪声信号n(m),m=1,2,...,M进行随机拟合,得到两道信号;然后根据独立成分分析算法的分离特性,运用独立成分分析算法对这两道信号进行分离处理,输出长度为M的峰度值不相等的两道信号;最后分别计算两道信号的峰度值,峰度值较大的信号记为x′(m),m=1,2,...,M,峰度值较小的信号记为n′(m),m=1,2,...,M;
步骤2:判断x′(m)与x(m)是否反相;如果反相,则置相位判别因子α=-1,否则置α=1;得到信号α·x′(m)(即若反相则进行自动反相校正,如果不反相,则保持不变);
步骤3:首先运用完全总体经验模态算法对α·x′(m)进行分解,得到P个信号分量,记为yp(m),p=1,...,P;m=1,...,M;然后计算各信号分量的峰度值,分别记为kp,p=1,...,P;再计算信号n′(m)的峰度值,记为k;最后以k作为阈值,将峰度值kp大于k的所有信号分量累加,作为原始回波信号x(m)去噪后的信号,记为
进一步地,所述步骤1中,零均值白噪声信号n(m)的幅度和方差任意设置。
进一步地,所述步骤1中,随机拟合的方法为:首先将x(m)和n(m)组合成2×M的矩阵,记为K,该矩阵的第一行为x(m),第二行为n(m),M表示x(m)的长度;然后生成一个2×2的矩阵记为L;最后将L左乘K,得到矩阵U,其第一行数据记为U1,第二行数据记为U2;U1和U2即为随机拟合得到的两道信号。
进一步地,所述步骤1和步骤3中,记某道信号为r(m),(m=1,2…,M),其峰度值计算公式如下:
其中,代表r(m),(m=1,2…,M)的均值。
进一步地,所述步骤2中,根据max(|x′(m)-x(m)|)>max(|x(m)|),m=0,...,M-1是否成立,判断x′(m)与x(m)是否反相;如果成立,说明x′(m)与x(m)反相;否则说明x′(m)与x(m)不反相;式中,|·|表示对一维数组x′(m)-x(m)中各元素取绝对值,max(|x(m)|)表示取一维数组|x(m)|中各元素的最大值。
有益效果:
本发明提出了一种基于自动反相校正和峰度值比较的GPR回波信号的去噪算法,针对独立成分分析算法分离后信号相位不定性问题,设计了相位判别因子,实现了独立成分分析算法分解后的信号相位自动判别和校正。针对完全总体经验模态分解算法分解后,各IMF分量需要人工判别的问题,设计了基于峰度值比较的IMF分量自动筛选,阈值选择为ICA算法分离出来的噪声信号的峰度值。本发明避免了了独立成分分析算法分解后的相位不定性,且在总体经验模态分解算法分解后无需传统的人工方式进行各IMF分量的筛选,提高了计算效率和去噪效果。
附图说明
图1示出了本方法的流程图。
图2示出了GPR正演模型图。
图3示出了图2正演模拟得到的GPR无噪回波信号图。
图4示出了人为加入噪声后模拟得到的含噪GPR信号图和产生的等长度的随机噪声信号;图4(a)为含噪GPR信号,图4(b)为根据图4(a)含噪GPR信号随机产生的等长度的随机噪声信号。
图5示出了图4所示的两个信号随机拟合后得到的两道信号,5(a)和(b)分别为随机拟合后得到的两道信号。
图6示出了图5所示的两道信号经过ICA算法后分离出的两道信号。
图7示出了图6中的第一道信号经过CEEMD分解后的各IMF分量波形图,为图7(a)为IMF5~IMF8,7(b)为IMF9~IMF12,7(c)为IMF13~IMF15。
图8示出了各IMF分量的峰度值计算结果和ICA分解得到的噪声信号的峰度值阈值。
图9示出了图7中的各IMF分量采用峰度值阈值比较后重构的信号。
图10示出了本方法和常规算法的去噪误差和信噪比的变化曲线对比图
具体实施方式
以下将结合附图和具体实施例对本发明做进一步详细说明。
GPR正演模拟如图2所示,地下介质有三层,第一二层的厚度均为15cm,第一二三层的相对介电常数分别为20,10和15,一个半径为4cm的金属管埋在第二层中,金属管的中心距第二层上表面6cm,发射和接收天线位于金属管中心的正上方,距地表的高度为5cm。发射信号为中心频率为900MHz的Ricker子波,时窗取为40ns,采用时域有限差分法仿真得到接收天线所接收的地下区域的散射回波,如图3所示,时长为40ns,采样点数为6784,该信号为无噪的回波信号s(m)。为模拟生成含噪的GPR信号,人为地加入了随机白噪声,白噪声的幅值和方差分别为4.3328和1.0513,加入噪声后的GPR回波信号如图4(a)所示,SNR=19。下面对该信号运用本发明的去噪算法进行去噪处理。根据图4(a)含噪GPR信号随机产生的等长度的随机噪声信号如图4(b)所示,再将图4(a)含噪GPR信号和图4(b)随机噪声信号进行随机拟合,即与一个2×2的矩阵相乘,得到两道混合信号,随机拟合结果如图5(a)和(b)所示。图5(a)和(b)输入ICA算法进行处理,得到两道输出信号,该两道信号的峰度值为81.6651和2.9763,峰度值高的输出信号记为x′(m),如图6(a)所示,峰度值低的输出信号记为n′(m),如图6(b)所示,对图6(a)信号进行判断是否信号反相,如果反相则对信号进行自动相位校正,如果未反相则信号不变。此处,经过公式判断,未发生反相,则α=1。再对α·x′(m)进行CEEMD算法的处理,得到图7所示的15个IMF分量。对图7中的每一个IMF分量计算峰度值,峰度值分别为[2.1453,2.1983,2.7546,2.8556,75.0391,41.8097,25.3151,16.4506,9.3599,7.3725,3.2763,3.7535,3.9701,1.6372,1.7456],如图8所示。如前所述,ICA算法分离出来的噪声信号n′(m),如图6(b)所示,其峰度值为2.9763。将图7中的每一个IMF分量的峰度值与图6(b)噪声信号的峰度值进行比较,将峰度值低于2.9763的IMF分量剔除,即剔除IMF1,IMF2,IMF3,IMF4,IMF14,IMF15分量,保留IMF5-IMF13分量。将IMF5-IMF13分量进行累加得到的重构信号作为去噪后的GPR回波信号z(m),如图9所示。将该信号与图3所示的无噪GPR回波信号相对比,并运用公式计算误差,得到均方误差为0.001085。
为定量评估该方法的去噪算法与常规CEEMD的去噪算法的性能,生成了不同SNR的含噪GPR回波信号,SNR的设置范围为[0,20],点数为21点。每一SNR情况下,根据图3的无噪的GPR回波,生成对应幅度的随机噪声信号,并将二者相加,作为含噪的GPR回波。分别运用本专利申请的算法和常规CEEMD算法对该回波进行去噪处理,得到各自的去噪回波。再分别与图3所示的无噪GPR回波进行对比分析,使用公式计算均方误差值。遍历整个SNR设置区间,得到各SNR情况下的均方误差曲线,如图10所示。
从图10中可见,本专利申请的去噪算法和常规CEEMD去噪算法随着SNR的增大,均方误差都随之降低。但在同一SNR情况下,本专利申请的去噪算法的均方误差更低,去噪效果更好。而且,本专利申请的去噪算法无需常规CEEMD去噪算法的人工方式剔除噪声分量的步骤,全程自动处理,效率更高。
Claims (2)
1.一种基于自动反相校正和峰度值比较的探地雷达回波信号去噪算法,其特征在于,包括以下步骤:
步骤1:首先将探地雷达单道含噪原始回波信号x(m),m=1,2,...,M与相同长度的零均值白噪声信号n(m),m=1,2,...,M进行随机拟合,得到两道信号;然后运用独立成分分析算法对这两道信号进行分离处理,输出长度为M的峰度值不相等的两道信号;最后分别计算两道信号的峰度值,峰度值较大的信号记为x′(m),m=1,2,...,M,峰度值较小的信号记为n′(m),m=1,2,...,M;
步骤2:判断x′(m)与x(m)是否反相;如果反相,则置相位判别因子α=-1,否则置α=1;得到信号α·x′(m);
步骤3:首先运用完全总体经验模态算法对α·x′(m)进行分解,得到P个信号分量,记为yp(m),p=1,...,P;m=1,...,M;然后计算各信号分量的峰度值,分别记为kp,p=1,...,P;再计算信号n′(m)的峰度值,记为k;最后以k作为阈值,将峰度值kp大于k的所有信号分量累加,作为原始回波信号x(m)去噪后的信号,记为
所述步骤1中,随机拟合的方法为:首先将x(m)和n(m)组合成2×M的矩阵,记为K,该矩阵的第一行为x(m),第二行为n(m),M表示x(m)的长度;然后生成一个2×2的矩阵记为L;最后将L左乘K,得到矩阵U,其第一行数据记为U1,第二行数据记为U2;U1和U2即为随机拟合得到的两道信号;
所述步骤1和步骤3中,记某道信号为r(m),m=1,2…,M,其峰度值计算公式如下:
其中,代表r(m),m=1,2…,M的均值;
所述步骤2中,根据max(|x′(m)-x(m)|)>max(|x(m)|),m=0,...,M-1是否成立,判断x′(m)与x(m)是否反相;如果成立,说明x′(m)与x(m)反相;否则说明x′(m)与x(m)不反相;式中,|·|表示对一维数组x′(m)-x(m)中各元素取绝对值,max(|x(m)|)表示取一维数组|x(m)|中各元素的最大值。
2.根据权利要求1所述的基于自动反相校正和峰度值比较的探地雷达回波信号去噪算法,其特征在于,所述步骤1中,零均值白噪声信号n(m)的幅度和方差任意设置。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201710296899.XA CN107121705B (zh) | 2017-04-28 | 2017-04-28 | 一种基于自动反相校正和峰度值比较的探地雷达回波信号去噪算法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201710296899.XA CN107121705B (zh) | 2017-04-28 | 2017-04-28 | 一种基于自动反相校正和峰度值比较的探地雷达回波信号去噪算法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN107121705A CN107121705A (zh) | 2017-09-01 |
CN107121705B true CN107121705B (zh) | 2018-10-12 |
Family
ID=59726035
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201710296899.XA Active CN107121705B (zh) | 2017-04-28 | 2017-04-28 | 一种基于自动反相校正和峰度值比较的探地雷达回波信号去噪算法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN107121705B (zh) |
Families Citing this family (3)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN106442635A (zh) * | 2016-09-22 | 2017-02-22 | 北京林业大学 | 一种基于雷达波的树木内部结构层位识别的方法 |
CN112285791B (zh) * | 2020-03-30 | 2022-01-25 | 中国科学院地质与地球物理研究所 | 一种探地雷达电磁场数据处理方法、装置以及介质 |
CN116299562B (zh) * | 2023-05-26 | 2023-08-04 | 中国海洋大学 | 一种高度计测距电离层误差校正滤波处理方法 |
Family Cites Families (6)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
DE602006016451D1 (de) * | 2005-07-28 | 2010-10-07 | Exxonmobil Upstream Res Co | Verfahren zur wavelet-entrauschung von elektromagnetischen vermessungsdaten gesteuerter quellen |
CN101822548A (zh) * | 2010-03-19 | 2010-09-08 | 哈尔滨工业大学(威海) | 基于相关性分析与经验模态分解的超声信号去噪方法 |
CN104020498B (zh) * | 2014-05-16 | 2016-10-19 | 桂林电子科技大学 | 基于瞬时幅频参数的探地雷达信号预处理方法及装置 |
CN104330793A (zh) * | 2014-10-15 | 2015-02-04 | 河海大学 | 一种基于Hilbert变换的探地雷达成像方法 |
CN104459809B (zh) * | 2014-10-30 | 2017-04-12 | 吉林大学 | 一种基于独立成分分析的全波核磁共振信号噪声滤除方法 |
CN105700020A (zh) * | 2016-03-23 | 2016-06-22 | 中国石油天然气集团公司 | 一种地震数据随机噪声压制方法及装置 |
-
2017
- 2017-04-28 CN CN201710296899.XA patent/CN107121705B/zh active Active
Also Published As
Publication number | Publication date |
---|---|
CN107121705A (zh) | 2017-09-01 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN102879766B (zh) | 一种高频地波雷达弱目标检测跟踪方法和装置 | |
Baili et al. | GPR signal de-noising by discrete wavelet transform | |
Wang et al. | Noise suppressing and direct wave arrivals removal in GPR data based on Shearlet transform | |
Leandro et al. | GPR trace analysis for coastal depositional environments of southern Brazil | |
Chen et al. | Adaptive ground clutter reduction in ground-penetrating radar data based on principal component analysis | |
CN109581516B (zh) | 曲波域统计量自适应阈值探地雷达数据去噪方法及*** | |
Bao et al. | GPR data noise attenuation on the curvelet transform | |
Zhang et al. | Advanced signal processing method for ground penetrating radar feature detection and enhancement | |
CN107121705B (zh) | 一种基于自动反相校正和峰度值比较的探地雷达回波信号去噪算法 | |
CN110133643B (zh) | 植物根系探测方法及装置 | |
CN108961181B (zh) | 一种基于shearlet变换的探地雷达图像去噪方法 | |
CN106199532A (zh) | 基于混合傅立叶‑小波分析的探地雷达信号降噪方法 | |
CN102495402A (zh) | 在钢筋强反射回波干扰下的道面灾害目标检测识别方法 | |
Hao et al. | An automated GPR signal denoising scheme based on mode decomposition and principal component analysis | |
CN113887398A (zh) | 一种基于变分模态分解和奇异谱分析的gpr信号去噪方法 | |
CN109085547B (zh) | 一种表层穿透雷达回波信号的去噪方法和相关装置 | |
CN108983158B (zh) | 一种基于Hankel矩阵奇异值分解的探地雷达噪声抑制方法 | |
Devapal et al. | Discontinuity adaptive SAR image despeckling using curvelet-based BM3D technique | |
Agrawal et al. | GPR data analysis of weak signals using modified S-transform | |
Kumar et al. | Removal of clutter and random noise for gpr images | |
CN112666552B (zh) | 一种探地雷达数据背景杂波自适应去除方法 | |
CN115113163A (zh) | 一种探地雷达多分辨率低秩稀疏分解杂波抑制方法 | |
CN101930605A (zh) | 基于二维混合变换的sar图像目标提取方法及*** | |
Xu et al. | Signal filtering based on wavelet transform and its application in ground penetrating radar | |
Ouadfeul et al. | Noise attenuation from GPR data using wavelet transform and artificial neural network |
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 |