CN105353373B - 一种基于Hough变换探地雷达目标提取方法和装置 - Google Patents
一种基于Hough变换探地雷达目标提取方法和装置 Download PDFInfo
- Publication number
- CN105353373B CN105353373B CN201510942842.3A CN201510942842A CN105353373B CN 105353373 B CN105353373 B CN 105353373B CN 201510942842 A CN201510942842 A CN 201510942842A CN 105353373 B CN105353373 B CN 105353373B
- Authority
- CN
- China
- Prior art keywords
- scan
- signal
- echo
- hyperbola
- wavelet
- 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
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
- G01S13/00—Systems using the reflection or reradiation of radio waves, e.g. radar systems; Analogous systems using reflection or reradiation of waves whose nature or wavelength is irrelevant or unspecified
- G01S13/88—Radar or analogous systems specially adapted for specific applications
- G01S13/885—Radar or analogous systems specially adapted for specific applications for ground probing
Landscapes
- Physics & Mathematics (AREA)
- Engineering & Computer Science (AREA)
- Radar, Positioning & Navigation (AREA)
- Remote Sensing (AREA)
- Electromagnetism (AREA)
- Computer Networks & Wireless Communication (AREA)
- General Physics & Mathematics (AREA)
- Radar Systems Or Details Thereof (AREA)
Abstract
本发明公开了一种基于Hough变换探地雷达目标提取方法和装置。该方法通过自主研发的探地雷达样机采集目标回波信号，并依据雷达回波特点进行零点对齐，然后采用滑动平移滤波器降低雷达样机电源噪声造成的回波图像噪声，在此基础上采用小波变换来分析回波的噪声特征，并通过设定阈值来检测目标回波的波峰和波谷等重要信息，最后采用Hough变换来提取目标的双曲线，从而获得目标的空间位置信息。本发明通过小波变换处理雷达回波信号，增强了信号中的有用信息，为后续处理的可靠性从根本上做了保障，同时改进的双曲线Hough变换有效提高了提取精度，极大地降低了计算量，适合工程应用。
Description
技术领域
本发明涉及探地雷达地下目标识别与提取技术领域，尤其涉及雷达信号处理以及数字图像处理技术领域，特别涉及一种基于Hough变换探地雷达目标提取方法和装置。
背景技术
探地雷达是一种新型的无损探测浅层地下目标的探测设备。它通过发射天线发射高频宽带电磁波，接收天线接收回波信号，最后通过对回波的分析推断地下介质的结构。由于地下不同介质的电磁性质差异，雷达波到达介质分界面时产生强烈的回波，根据回波信号的旅程时间、幅度、相位信息，推断介质的结构及性质等。
雷达脉冲信号在发射到接收的过程中，主要的信号衰减是来自于地层的两次传播、分界面的反射、收发端的能量损耗。由于收发端的能量损耗是固定不变且一定范围内可控的，故可以对地层中雷达信号的传输、散射和反射进行分析。
探地雷达在一个水平表面上工作。设x为雷达的水平位置，电磁波在介质中传播，遇到目标反射回来由接收天线接收，这之间的旅程时间设为t，当雷达处于目标上方时，设x＝x0,t＝t0。那么推导出的x-t方程为：
针对雷达回波信号列构成的二维灰度图像，检测出上式所示的双曲线，并确定双曲线的参数便是最亟需解决的问题。
发明内容
本发明的目的在于提供一种基于Hough变换探地雷达目标提取方法和装置，以解决探地雷达目标提取与识别的相关问题，通过信号处理、数字图像处理、图形识别等一系列方法来获得地下目标的空间位置信息。
一种基于Hough变换探地雷达目标提取方法，获取N列A-Scan回波信号构成二维灰度图像B-Scan，再通过分析处理B-Scan的数据来获取目标位置信息，具体包括以下步骤：
步骤S1、对获取的N列A-Scan回波信号数据预处理与归一化；
步骤S2、滑动平移滤波；
步骤S3、模极大值判别波峰波谷；
步骤S4、基于Hough变换提取双曲线；
步骤S5、实测物体位置计算。
其中，所述步骤S1、对获取的N列A-Scan回波信号数据预处理与归一化；具体为：
将A-Scan信号数据中前一部分无效数据点去除，搜索待处理数据的最大值
与最小值，并将数据归一化到[0,1]之内以备后续处理，归一化公式为：
X″＝(X-min(X))/(max(X)-min(X))
其中，X″代表归一化后的回波数据值，X代表归一化前的回波数据值，max(X)代表回波数据中的最大值，min(X)代表回波数据中的最小值。
其中，所述步骤S2、滑动平移滤波；具体为：
对于步骤S1中得到的A-Scan数据，采用16组相邻A-Scan信号取滑动平均的方法进行平滑处理：
X(:，i)＝Σ0≤k＜16X(:，i+k)/16，0≤i＜N-16
其中，N为二维图像B-Scan中回波信号的总列数。
其中，所述步骤S3、模极大值判别波峰波谷；具体为：
(a)选定coif3小波进行小波变换；
(b)将各列A-Scan原始信号进行coif3小波分解，对小波系数进行处理之后重构信号；
(c)通过设定高阈值和低阈值来确定峰值点出现的区域，并使用导数过零点的方法找到峰值点准确位置；
(d)将提取到的二维波峰点二值化为灰度图，其中，峰值点为1，其他为0。
其中，所述步骤S4、基于Hough变换提取双曲线；具体为：
针对双曲线方程：
(1)匹配A，B；
根据双曲线存在一对渐近线：
可以得到：
其中(x1,y1),(x2,y2)为图像空间的两个点，实际有效点为双曲线两支所在位置，又有则得到A,B之间的关系；
(2)匹配a，b；
根据双曲线顶点曲线类似圆弧的原理提取实轴与虚轴长度，具体方法：抽取某个区域三个不同块的点，移动该区域获得不同组抽样值，计算三个点组成的弧的圆心位置，用Hough变换搜索最多的圆心点的位置并依据该位置来确定双曲线的参数；
设求得的圆方程为：
(x-x1)2+(y-y1)2＝r2
圆的参数与双曲线的参数关系为：
a＝x1,b＝y1-r-B,r＝B2/A；
(3)求解双曲线参数；
联立得以下方程组并求解即可得出双曲线参数：
其中，所述步骤S5、实测物体位置计算；具体为：
探地雷达目标构成的双曲线方程为：
采用步骤S4的方法确定x0，t0以及v，则目标在坐标系中的位置为(x0,vt0)。
其中，所述(b)将各列A-Scan原始信号进行coif3小波分解，对小波系数进行处理之后重构信号；具体为：
1)选定coif3小波进行小波变换，对每一列A-Scan原始信号s进行5次coif3小波分解，得到第五次分解的低频信息α5，以及相应每层分解的高频细节信息d1，d2，d3，d4，d5；
2)将能够反映信号大尺度上的交流分量的α5和能够反映信号的高频噪声信息的d1，d2的小波系数置为零；
3)保留d4的小波系数，同时将d3,d5的小波系数按一定比例缩小；
4)将剩余的小波系数进行逆小波变换，对原始信号进行重构。
一种基于Hough变换探地雷达目标提取装置，该装置用于获取N列A-Scan回波信号构成二维灰度图像B-Scan，再通过分析处理B-Scan的数据来获取目标位置信息，具体包括以下单元：
数据预处理单元，用于对获取的N列A-Scan回波信号数据预处理与归一化；
滤波除噪单元，用于滑动平移滤波；
波峰波谷获取单元，用于模极大值判别波峰波谷；
双曲线提取单元，用于基于Hough变换提取双曲线；
物体位置计算单元，用于实测物体位置计算。
其中，所述数据预处理单元，具体用于将A-Scan信号数据中前一部分无效数据点去除，搜索待处理数据的最大值与最小值，并将数据归一化到[0,1]之内以备后续处理，归一化公式为：
X″＝(X-min(X))/(max(X)-min(X))
其中，X″代表归一化后的回波数据值，X代表归一化前的回波数据值，max(X)代表回波数据中的最大值，min(X)代表回波数据中的最小值；
所述滤波除噪单元，具体用于对于数据预处理单元中得到的A-Scan数据，采用16组相邻A-Scan信号取滑动平均的方法进行平滑处理：
X(:，i)＝Σ0≤k＜16X(:，i+k)/16，0≤i＜N-16
其中，N为二维图像B-Scan中回波信号的总列数；
所述波峰波谷获取单元，具体用于选定coif3小波进行小波变换；将各列A-Scan原始信号进行coif3小波分解，对小波系数进行处理之后重构信号；通过设定高阈值和低阈值来确定峰值点出现的区域，并使用导数过零点的方法找到峰值点准确位置；将提取到的二维波峰点二值化为灰度图，其中，峰值点为1，其他为0；
所述双曲线提取单元，具体用于针对双曲线方程：
(1)匹配A，B；
根据双曲线存在一对渐近线：
可以得到：
其中(x1,y1),(x2,y2)为图像空间的两个点，实际有效点为双曲线两支所在位置，又有则得到A,B之间的关系；
(2)匹配a，b；
根据双曲线顶点曲线类似圆弧的原理提取实轴与虚轴长度，具体方法：抽取某个区域三个不同块的点，移动该区域获得不同组抽样值，计算三个点组成的弧的圆心位置，用Hough变换搜索最多的圆心点的位置并依据该位置来确定双曲线的参数；
设求得的圆方程为：
(x-x1)2+(y-y1)2＝r2
圆的参数与双曲线的参数关系为：
a＝x1,b＝y1-r-B,r＝B2/A；
(3)求解双曲线参数；
联立得以下方程组并求解即可得出双曲线参数：
所述物体位置计算单元，具体用于探地雷达目标构成的双曲线方程为：
采用双曲线提取单元的计算方法确定x0，t0以及v，则目标在坐标系中的位置为(x0,vt0)。
其中，所述将各列A-Scan原始信号进行coif3小波分解，对小波系数进行处理之后重构信号；具体为：
1)选定coif3小波进行小波变换，对每一列A-Scan原始信号s进行5次coif3小波分解，得到第五次分解的低频信息α5，以及相应每层分解的高频细节信息d1，d2，d3，d4，d5；
2)将能够反映信号大尺度上的交流分量的α5和能够反映信号的高频噪声信息的d1，d2的小波系数置为零；
3)保留d4的小波系数，同时将d3,d5的小波系数按一定比例缩小；
4)将剩余的小波系数进行逆小波变换，对原始信号进行重构。
有益效果：
本发明公开了一种基于Hough变换探地雷达目标提取方法和装置，所述方法通过获取N列A-Scan回波信号构成二维灰度图像B-Scan，再通过分析处理B-Scan的数据来获取目标位置信息，具体包括以下步骤：步骤S1、对获取的N列A-Scan回波信号数据预处理与归一化；步骤S2、滑动平移滤波；步骤S3、模极大值判别波峰波谷；步骤S4、基于Hough变换提取双曲线；步骤S5、实测物体位置计算。本发明解决了探地雷达目标提取与识别的相关问题，通过信号处理、数字图像处理、图形识别等一系列方法来获得地下目标的空间位置信息。本发明通过小波变换处理雷达回波信号，增强了信号中的有用信息，为后续处理的可靠性从根本上做了保障，同时改进的双曲线Hough变换有效提高了提取精度，极大地降低了计算量，适合工程应用。
附图说明
图1：本发明提供的一种基于Hough变换探地雷达目标提取方法的流程图；
图2：A-Scan数据叠加图；
图3：单行A-Scan图；
图4：原始B-Scan图；
图5：滑动平均B-Scan图；
图6：A-Scan信号小波分解图；
图7：原始信号与重构信号对比图；
图8：波峰二值化图；
图9：双曲线提取图；
图10：本发明提供的一种基于Hough变换探地雷达目标提取装置的结构图。
具体实施方式
下面结合附图和实施例对本发明技术方案作进一步详细具体的说明。
实施例1：
图1是本发明提供的一种基于Hough变换探地雷达目标提取方法的流程图。如图1所示，本发明所述的一种基于Hough变换探地雷达目标提取方法，获取N列A-Scan回波信号构成二维灰度图像B-Scan，再通过分析处理B-Scan的数据来获取目标位置信息，具体包括以下步骤：
步骤S1、对获取的N列A-Scan回波信号数据预处理与归一化；
步骤S2、滑动平移滤波；
步骤S3、模极大值判别波峰波谷；
步骤S4、基于Hough变换提取双曲线；
步骤S5、实测物体位置计算。
本发明解决了探地雷达目标提取与识别的相关问题，通过信号处理、数字图像处理、图形识别等一系列方法来获得地下目标的空间位置信息。
图2是A-Scan数据叠加图。如图2所示为B-Scan中每一列雷达回波信号A-Scan在同一坐标系内的叠加显示，本发明要处理的便是这部分数据。图2中横坐标显示A-Scan中点的标号，纵坐标显示为16位AD的采样值。每列图像有512个点，大概在40个点以后数据有明显的波形起伏，也就是有用的回波信号；前面的数据是无用的，首先将回波信号前面的无用数据点去除，由于matlab灰度图的显示范围是在[0,1]以内，因此应确定回波信号的最大值与最小值，并以此为基础将数据归一化到[0,1]之间，以备之后的图像显示和数据处理。因此去除前面30个无效数据点，剩余数据分布在[-30000,30000]之间。将这些数据归一化到[0,1]之间，以便之后的显示和数据处理。
所述步骤S1、对获取的N列A-Scan回波信号数据预处理与归一化；具体为：
将A-Scan信号数据中前一部分无效数据点去除，搜索待处理数据的最大值
与最小值，并将数据归一化到[0,1]之内以备后续处理，归一化公式为：
X″＝(X-min(x))/(max(X)-min(X))
其中，X″代表归一化后的回波数据值，X代表归一化前的回波数据值，max(X)代表回波数据中的最大值，min(X)代表回波数据中的最小值。
与所述步骤S1相同，这里同样是对每一列回波信号进行处理。由于雷达回波信号很弱，在通过程控放大和低通滤波之后仍然存在较多噪声，这直接导致合成的二维灰度图同样存在过多噪声。所以针对图3所示单列A-Scan回波信号进行滑动平移滤波。参照该雷达系统的参数，该雷达系统的脉冲扫描周期是256次/秒，同时雷达天线的拖动速度为0.1m/s，相邻两列信号的水平距离为：
0.1÷256＝0.00039m＝0.39mm
这个距离远小于雷达系统的水平分辨率。对16组相邻A-Scan信号，最大的水平距离为0.39mm×16＝6.24mm，对于实际运用的分辨率也在可接受范围内。因此，我们采用相邻16列归一化A-Scan信号取滑动平均的方法进行处理。
所述步骤S2、滑动平移滤波；具体为：
对于步骤S1中得到的A-Scan数据，采用16组相邻A-Scan信号取滑动平均的方法进行平滑处理：
X(:，i)＝Σ0≤k＜16X(:，i+k)/16，0≤i＜N-16
其中，N为二维图像B-Scan中回波信号的总列数。
最后得到的数据灰度图如图5所示，相比图4所示原始B-Scan灰度图，可以明显看到噪声的减小。
所述步骤S3、模极大值判别波峰波谷；具体为：
(a)选定coif3小波进行小波变换；
首先，小波基选取。针对实际情况，信号中会夹杂有高频噪声，造成信号小波变换的导数出现很多过零点；信号也夹杂了低频噪声，使得直接用阈值对信号进行判断很难准确地确定所有有效的波峰和波谷的位置，因此要对雷达回波信号进行小波变换。相同的信号选取不同的小波基进行小波变换，得到的小波有很大的不同。雷达系统的发射天线相当于一个高Q值的发射模块，信号大部分以某些频段的正弦分量的单脉冲形式通过天线耦合出去。信号进入地层，遇到地下分界面发生反射，最后返回到接收装置转换为电信号。这个过程相当于信号经过一个线性时不变系统后的结果，而该传输过程可以看成是经过一个滤波器，得到的回波信号基本成分不变，由很多正弦分量组成，只是幅值、相位有些差异。雷达发射天线和雷达信号传输特性使得回波信号的能量在频域比较集中，我们希望选取的母小波也是频域能量集中的小波。
Coif3小波在时域上有很好的单峰性、对称性、类正弦性。在频域上，该小波的能量也非常集中，能使信号的小波分解系数相对集中。同时小波在单峰下降后没有出现振铃现象，也没有出现突变现象，能够很好的在不同尺度上反应回波信号的不同分量。因此，选用coif3小波对信号进行小波变换，来完成波峰、波谷的检测。
(b)将各列A-Scan原始信号进行coif3小波分解，对小波系数进行处理之后重构信号；具体为：
1)选定coif3小波进行小波变换，对每一列A-Scan原始信号s进行5次coif3小波分解，得到第五次分解的低频信息α5，以及相应每层分解的高频细节信息d1，d2，d3，d4，d5；
先对每一列回波信号进行5次coif3小波分解，结果如图6所示。其中，s为原始信号。α5是第五次分解的低频信息，d5,d4,d3,d2,d1是相应的每层分解的高频细节信息。其中，d4信号最能真实的反应原始信号s中波峰和波谷信息；α5反应了信号的大尺度上的交流分量；d1,d2反映了信号的高频噪声信息。
2)将能够反映信号大尺度上的交流分量的α5和能够反映信号的高频噪声信息的d1，d2的小波系数置为零；
3)保留d4的小波系数，同时将d3,d5的小波系数按一定比例缩小；
4)将剩余的小波系数进行逆小波变换，对原始信号进行重构。重构信号如图7所示。
(c)通过设定高阈值和低阈值来确定峰值点出现的区域，并使用导数过零点的方法找到峰值点准确位置；
(d)将提取到的二维波峰点二值化为灰度图，其中，峰值点为1，其他为0。将重构信号A-Scan合成B-Scan灰度图之后，进行二维波峰点二值化结果如图8所示。
所述步骤S4、基于Hough变换提取双曲线；根据直线Hough变换类比推理，扩展到广义Hough变换的一种—双曲线Hough变换。对于双曲线方程：
根据目标双曲线的尺度大小，对图像空间进行分块。块内不进行抽样，只对块与块之间抽取样点进行参数空间计算。这样可以极大地减少图像空间的抽样组数，同时也不会丢失有效的信息。
对于匹配参数空间4个参数(a,b,A,B)，转化为先匹配可能存在的所有两个参数(a,b)或(A,B)的运算：
(1)匹配A，B；
根据双曲线存在一对渐近线：
可以得到：
其中(x1,y1),(x2,y2)为图像空间的两个点，实际有效点为双曲线两支所在位置，又有则得到A,B之间的关系；
(2)匹配a，b；
根据双曲线顶点曲线类似圆弧的原理提取实轴与虚轴长度，具体方法：抽取某个区域三个不同块的点，移动该区域获得不同组抽样值，计算三个点组成的弧的圆心位置，用Hough变换搜索最多的圆心点的位置作为双曲线的交点位置，并依据该位置来确定双曲线的参数；
设求得的圆方程为：
(x-x1)2+(y-y1)2＝r2
圆的参数与双曲线的参数关系为：
a＝x1,b＝y1-r-B,r＝B2/A；
(3)求解双曲线参数；
联立得以下方程组并求解即可得出双曲线参数：
所述步骤S5、实测物体位置计算；具体为：
探地雷达目标构成的双曲线方程为：
采用步骤S4的方法确定x0，t0以及v，则目标在坐标系中的位置为(x0,vt0)。
根据图8所示波峰二值化图可知，双曲线的水平跨度是2000，垂直跨度是200，双曲线的渐近线大概在水平30度左右，左右对称没有明显的旋转误差。因此我们选取一些初设参数进行分析：
分块处理：在水平方向以100为跨度，垂直方向以20为跨度进行分块。全幅图像分层大概25×30块。
|k|的范围选择：由于图中很多水平直线，同时渐近线的斜率在水平30度左右，故选取渐进线的倾斜角在10°～45°之间，因此|k|应该满足：
0.17＜|k|＜1，经过Hough变换后确定|k|＝0.18。
顶点弧抽样块范围选择：选取一个2×6的模板进行滑动抽样，抽取模板内的三个点进行顶点曲率圆参数计算。
最终得出最明显的双曲线如图9所示，其表达式为：
由于雷达回波和起始采样不对齐等问题，使得每一列回波信号的第一个点并不是回波信号的起始点。故我们搜索并将直达波的位置作为起始位置。直达波可以看作是发射脉冲直接通过天线外壳反射回来的信号。当天线紧贴地面时，直达波旅程时间可以看作脉冲在空气中传播的时间，故用回波信号时间减去直达波旅程时间等价为脉冲信号在土层中经过反射往返传播所用的时间。根据前面的滑动平均滤波的二维回波图，可以确定直达波在第68个点的位置。令x-1878＝0,得到双曲线顶点纵坐标为：
y＝108+30＝138
回波双曲线顶点位置在138个点的位置。故Δn＝138-68＝70。然而每列雷达回波信号相邻两个点的等效时间差为0.1ns，雷达信号在介电常数为εr介质中的传播速度为：
其中c为真空光速为0.3m/ns，根据查表可知肥沃湿润土壤的介电常数εr约为6～10。故实际物理埋藏深度为：
由于实际土壤介电常数εr没有实际做实验测量，计算的结果0.35m和实际深度0.4m相比是处于误差范围内的，这就验证了本发明方法的有效性。
本发明通过小波变换处理雷达回波信号，增强了信号中的有用信息，为后续处理的可靠性从根本上做了保障，同时改进的双曲线Hough变换有效提高了提取精度，极大地降低了计算量，适合工程应用。
实施例2：
本发明所述的装置实施例2与方法实施例1属于统一技术构思，在装置实施例2中未详尽描述的内容，请参见方法实施例1。
图10：本发明提供的一种基于Hough变换探地雷达目标提取装置的结构图。如图10所示，本发明所述的一种基于Hough变换探地雷达目标提取装置，该装置用于获取N列A-Scan回波信号构成二维灰度图像B-Scan，再通过分析处理B-Scan的数据来获取目标位置信息，具体包括以下单元：
数据预处理单元，用于对获取的N列A-Scan回波信号数据预处理与归一化；
滤波除噪单元，用于滑动平移滤波；
波峰波谷获取单元，用于模极大值判别波峰波谷；
双曲线提取单元，用于基于Hough变换提取双曲线；
物体位置计算单元，用于实测物体位置计算。
所述数据预处理单元，具体用于将A-Scan信号数据中前一部分无效数据点去除，搜索待处理数据的最大值与最小值，并将数据归一化到[0,1]之内以备后续处理，归一化公式为：
X″＝(X-min(X))/(max(X)-min(X))
其中，X″代表归一化后的回波数据值，X代表归一化前的回波数据值，max(X)代表回波数据中的最大值，min(X)代表回波数据中的最小值。
图2是A-Scan数据叠加图。如图2所示为B-Scan中每一列雷达回波信号A-Scan在同一坐标系内的叠加显示，本发明要处理的便是这部分数据。图2中横坐标显示A-Scan中点的标号，纵坐标显示为16位AD的采样值。每列图像有512个点，大概在40个点以后数据有明显的波形起伏，也就是有用的回波信号；前面的数据是无用的，首先将回波信号前面的无用数据点去除，由于matlab灰度图的显示范围是在[0,1]以内，因此应确定回波信号的最大值与最小值，并以此为基础将数据归一化到[0,1]之间，以备之后的图像显示和数据处理。因此去除前面30个无效数据点，剩余数据分布在[-30000,30000]之间。将这些数据归一化到[0,1]之间，以便之后的显示和数据处理。
所述滤波除噪单元，具体用于对于数据预处理单元中得到的A-Scan数据，采用16组相邻A-Scan信号取滑动平均的方法进行平滑处理：
X(:，i)＝Σ0≤k＜16X(:，i+k)/16，0≤i＜N-16
其中，N为二维图像B-Scan中回波信号的总列数。与发明内容所述数据预处理单元相同，这里同样是对每一列回波信号进行处理。由于雷达回波信号很弱，在通过程控放大和低通滤波之后仍然存在较多噪声，这直接导致合成的二维灰度图同样存在过多噪声。所以针对图3所示单列A-Scan回波信号进行滑动平移滤波。参照该雷达系统的参数，该雷达系统的脉冲扫描周期是256次/秒，同时雷达天线的拖动速度为0.1m/s，相邻两列信号的水平距离为：
0.1m/s÷256/s＝0.00039m＝0.39mm
这个距离远小于雷达系统的水平分辨率。对16组相邻A-Scan信号，最大的水平距离为0.39mm×16＝6.24mm，对于实际运用的分辨率也在可接受范围内。因此，我们采用相邻16列归一化A-Scan信号取滑动平均的方法进行处理。
最后得到的数据灰度图如图5所示，相比图4所示原始B-Scan灰度图，可以明显看到噪声的减小。
所述波峰波谷获取单元，具体用于(1)选定coif3小波进行小波变换；(2)将各列A-Scan原始信号进行coif3小波分解，对小波系数进行处理之后重构信号；(3)通过设定高阈值和低阈值来确定峰值点出现的区域，并使用导数过零点的方法找到峰值点准确位置；(4)将提取到的二维波峰点二值化为灰度图，其中，峰值点为1，其他为0。将重构信号A-Scan合成B-Scan灰度图之后，进行二维波峰点二值化结果如图8所示。
首先，小波基选取。针对实际情况，信号中会夹杂有高频噪声，造成信号小波变换的导数出现很多过零点；信号也夹杂了低频噪声，使得直接用阈值对信号进行判断很难准确地确定所有有效的波峰和波谷的位置，因此要对雷达回波信号进行小波变换。相同的信号选取不同的小波基进行小波变换，得到的小波有很大的不同。雷达系统的发射天线相当于一个高Q值的发射模块，信号大部分以某些频段的正弦分量的单脉冲形式通过天线耦合出去。信号进入地层，遇到地下分界面发生反射，最后返回到接收装置转换为电信号。这个过程相当于信号经过一个线性时不变系统后的结果，而该传输过程可以看成是经过一个滤波器，得到的回波信号基本成分不变，由很多正弦分量组成，只是幅值、相位有些差异。雷达发射天线和雷达信号传输特性使得回波信号的能量在频域比较集中，我们希望选取的母小波也是频域能量集中的小波。
Coif3小波在时域上有很好的单峰性、对称性、类正弦性。在频域上，该小波的能量也非常集中，能使信号的小波分解系数相对集中。同时小波在单峰下降后没有出现振铃现象，也没有出现突变现象，能够很好的在不同尺度上反应回波信号的不同分量。因此，选用coif3小波对信号进行小波变换，来完成波峰、波谷的检测。
所述将各列A-Scan原始信号进行coif3小波分解，对小波系数进行处理之后重构信号；具体为：
1)选定coif3小波进行小波变换，对每一列A-Scan原始信号s进行5次coif3小波分解，得到第五次分解的低频信息α5，以及相应每层分解的高频细节信息d1，d2，d3，d4，d5；
先对每一列回波信号进行5次coif3小波分解，结果如图6所示。其中，s为原始信号。α5是第五次分解的低频信息，d5,d4,d3,d2,d1是相应的每层分解的高频细节信息。其中，d4信号最能真实的反应原始信号s中波峰和波谷信息；α5反应了信号的大尺度上的交流分量；d1,d2反映了信号的高频噪声信息。
2)将能够反映信号大尺度上的交流分量的α5和能够反映信号的高频噪声信息的d1，d2的小波系数置为零；
3)保留d4的小波系数，同时将d3,d5的小波系数按一定比例缩小；
4)将剩余的小波系数进行逆小波变换，对原始信号进行重构。重构信号如图7所示。
所述双曲线提取单元，具体用于针对双曲线方程：
根据目标双曲线的尺度大小，对图像空间进行分块。块内不进行抽样，只对块与块之间抽取样点进行参数空间计算。这样可以极大地减少图像空间的抽样组数，同时也不会丢失有效的信息。
对于匹配参数空间4个参数(a,b,A,B)，转化为先匹配可能存在的所有两个参数(a,b)或(A,B)的运算：
(1)匹配A，B；
根据双曲线存在一对渐近线：
可以得到：
其中(x1,y1),(x2,y2)为图像空间的两个点，实际有效点为双曲线两支所在位置，又有则得到A,B之间的关系；
(2)匹配a，b；
根据双曲线顶点曲线类似圆弧的原理提取实轴与虚轴长度，具体方法：抽取某个区域三个不同块的点，移动该区域获得不同组抽样值，计算三个点组成的弧的圆心位置，用Hough变换搜索最多的圆心点的位置并依据该位置来确定双曲线的参数；
设求得的圆方程为：
(x-x1)2+(y-y1)2＝r2
圆的参数与双曲线的参数关系为：
a＝x1,b＝y1-r-B,r＝B2/A；
(3)求解双曲线参数；
联立得以下方程组并求解即可得出双曲线参数：
所述物体位置计算单元，具体用于探地雷达目标构成的双曲线方程为：
采用双曲线提取单元的计算方法确定x0，t0以及v，则目标在坐标系中的位置为(x0,vt0)。
以上所述实施例及应用场景仅为本发明的较佳实施例及应用场景而已，并不用以限制本发明，凡在本发明的精神和原则之内所作的任何修改、等同替换和改进及其在其他领域及场景的应用，均应包含在本发明的保护范围之内。
Claims (8)
1.一种基于Hough变换探地雷达目标提取方法，其特征在于，获取N列A-Scan回波信号构成二维灰度图像B-Scan，再通过分析处理B-Scan的数据来获取目标位置信息，具体包括以下步骤：
步骤S1、对获取的N列A-Scan回波信号数据预处理与归一化；
步骤S2、滑动平移滤波；
步骤S3、模极大值判别波峰波谷；
步骤S4、基于Hough变换提取双曲线；
步骤S5、实测物体位置计算；
所述步骤S3、模极大值判别波峰波谷；具体为：
(a)选定coif3小波进行小波变换；
(b)将各列A-Scan原始信号进行coif3小波分解，对小波系数进行处理之后重构信号；
(c)通过设定高阈值和低阈值来确定峰值点出现的区域，并使用导数过零点的方法找到峰值点准确位置；
(d)将提取到的二维波峰点二值化为灰度图，其中，峰值点为1，其他为0，
所述步骤S4、基于Hough变换提取双曲线；具体为：
针对双曲线方程：
(1)匹配A，B；
根据双曲线存在一对渐近线：
可以得到：
其中(x1,y1),(x2,y2)为图像空间的两个点，实际有效点为双曲线两支所在位置，又有则得到A,B之间的关系；
(2)匹配a，b；
根据双曲线顶点曲线类似圆弧的原理提取实轴与虚轴长度，具体方法：抽取某个区域三个不同块的点，移动该区域获得不同组抽样值，计算三个点组成的弧的圆心位置，用Hough变换搜索最多的圆心点的位置并依据该位置来确定双曲线的参数；
设求得的圆方程为：
(x-x1)2+(y-y1)2＝r2
圆的参数与双曲线的参数关系为：
a＝x1,b＝y1-r-B,r＝B2/A；
(3)求解双曲线参数；
联立得以下方程组并求解即可得出双曲线参数：
2.根据权利要求1所述的一种基于Hough变换探地雷达目标提取方法，其特征在于，所述步骤S1、对获取的N列A-Scan回波信号数据预处理与归一化；具体为：
将A-Scan信号数据中前一部分无效数据点去除，搜索待处理数据的最大值与最小值，并将数据归一化到[0,1]之内以备后续处理，归一化公式为：
X″＝(X-min(x))/(max(X)-min(X))
其中，X″代表归一化后的回波数据值，X代表归一化前的回波数据值，max(X)代表回波数据中的最大值，min(X)代表回波数据中的最小值。
3.根据权利要求2所述的一种基于Hough变换探地雷达目标提取方法，其特征在于，所述步骤S2、滑动平移滤波；具体为：
对于步骤S1中得到的A-Scan数据，采用16组相邻A-Scan信号取滑动平均的方法进行平滑处理：
X(：，i)＝∑0≤k＜16X(：，i+k)/16，0≤i＜N-16
其中，N为二维图像B-Scan中回波信号的总列数。
4.根据权利要求1所述的一种基于Hough变换探地雷达目标提取方法，其特征在于，所述步骤S5、实测物体位置计算；具体为：
探地雷达目标构成的双曲线方程为：
采用步骤S4的方法确定x0，t0以及v，则目标在坐标系中的位置为(x0,vt0)。
5.根据权利要求1所述的一种基于Hough变换探地雷达目标提取方法，其特征在于，所述(b)将各列A-Scan原始信号进行coif3小波分解，对小波系数进行处理之后重构信号；具体为：
1)选定coif3小波进行小波变换，对每一列A-Scan原始信号s进行5次coif3小波分解，得到第五次分解的低频信息α5，以及相应每层分解的高频细节 信息d1，d2，d3，d4，d5；
2)将能够反映信号大尺度上的交流分量的α5和能够反映信号的高频噪声信息的d1，d2的小波系数置为零；
3)保留d4的小波系数，同时将d3,d5的小波系数按一定比例缩小；
4)将剩余的小波系数进行逆小波变换，对原始信号进行重构。
6.一种基于Hough变换探地雷达目标提取装置，其特征在于，该装置用于获取N列A-Scan回波信号构成二维灰度图像B-Scan，再通过分析处理B-Scan的数据来获取目标位置信息，具体包括以下单元：
数据预处理单元，用于对获取的N列A-Scan回波信号数据预处理与归一化；
滤波除噪单元，用于滑动平移滤波；
波峰波谷获取单元，用于模极大值判别波峰波谷；
双曲线提取单元，用于基于Hough变换提取双曲线；
物体位置计算单元，用于实测物体位置计算；
所述波峰波谷获取单元，具体用于选定coif3小波进行小波变换；将各列A-Scan原始信号进行coif3小波分解，对小波系数进行处理之后重构信号；通过设定高阈值和低阈值来确定峰值点出现的区域，并使用导数过零点的方法找到峰值点准确位置；将提取到的二维波峰点二值化为灰度图，其中，峰值点为1，其他为0，
所述双曲线提取单元，具体用于针对双曲线方程：
(1)匹配A，B；
根据双曲线存在一对渐近线：
可以得到：
其中(x1,y1),(x2,y2)为图像空间的两个点，实际有效点为双曲线两支所在位置，又有则得到A,B之间的关系；
(2)匹配a，b；
根据双曲线顶点曲线类似圆弧的原理提取实轴与虚轴长度，具体方法：抽取某个区域三个不同块的点，移动该区域获得不同组抽样值，计算三个点组成的弧的圆心位置，用Hough变换搜索最多的圆心点的位置并依据该位置来确定双曲线的参数；
设求得的圆方程为：
(x-x1)2+(y-y1)2＝r2
圆的参数与双曲线的参数关系为：
a＝x1,b＝y1-r-B,r＝B2/A；
(3)求解双曲线参数；
联立得以下方程组并求解即可得出双曲线参数：
7.根据权利要求6所述的一种基于Hough变换探地雷达目标提取装置，其特征在于，所述数据预处理单元，具体用于将A-Scan信号数据中前一部分无效数据点去除，搜索待处理数据的最大值与最小值，并将数据归一化到[0,1]之内 以备后续处理，归一化公式为：
X″＝(X-min(x))/(max(X)-min(X))
其中，X″代表归一化后的回波数据值，X代表归一化前的回波数据值，max(X)代表回波数据中的最大值，min(X)代表回波数据中的最小值；
所述滤波除噪单元，具体用于对于数据预处理单元中得到的A-Scan数据，采用16组相邻A-Scan信号取滑动平均的方法进行平滑处理：
X(：，i)＝∑0≤k＜16X(：，i+k)/16，0≤i＜N-16
其中，N为二维图像B-Scan中回波信号的总列数；
所述物体位置计算单元，具体用于探地雷达目标构成的双曲线方程为：
采用双曲线提取单元的计算方法确定x0，t0以及v，则目标在坐标系中的位置为(x0,vt0)。
8.根据权利要求6所述的一种基于Hough变换探地雷达目标提取装置，其特征在于，所述将各列A-Scan原始信号进行coif3小波分解，对小波系数进行处理之后重构信号；具体为：
1)选定coif3小波进行小波变换，对每一列A-Scan原始信号s进行5次coif3小波分解，得到第五次分解的低频信息α5，以及相应每层分解的高频细节信息d1，d2，d3，d4，d5；
2)将能够反映信号大尺度上的交流分量的α5和能够反映信号的高频噪声信息的d1，d2的小波系数置为零；
3)保留d4的小波系数，同时将d3,d5的小波系数按一定比例缩小；
4)将剩余的小波系数进行逆小波变换，对原始信号进行重构。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201510942842.3A CN105353373B (zh) | 2015-12-16 | 2015-12-16 | 一种基于Hough变换探地雷达目标提取方法和装置 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201510942842.3A CN105353373B (zh) | 2015-12-16 | 2015-12-16 | 一种基于Hough变换探地雷达目标提取方法和装置 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN105353373A CN105353373A (zh) | 2016-02-24 |
| CN105353373B true CN105353373B (zh) | 2018-04-20 |
Family
ID=55329374
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201510942842.3A Active CN105353373B (zh) | 2015-12-16 | 2015-12-16 | 一种基于Hough变换探地雷达目标提取方法和装置 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN105353373B (zh) |
Families Citing this family (11)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN107346023A (zh) * | 2017-07-04 | 2017-11-14 | 山东工商学院 | 一种基于矩特征的探地雷达双曲波显著性映射方法 |
| CN107358617B (zh) * | 2017-08-07 | 2020-03-20 | 中南大学 | 一种探地雷达记录剖面中目标顶点的检测方法 |
| CN108710115B (zh) * | 2018-05-14 | 2021-07-20 | 桂林电子科技大学 | 基于小波分解的探地雷达目标的检测方法 |
| CN108469609B (zh) * | 2018-06-11 | 2022-02-18 | 成都纳雷科技有限公司 | 一种用于雷达目标跟踪的检测信息滤波方法 |
| CN109283524B (zh) * | 2018-08-01 | 2020-11-06 | 西安交通大学 | 一种用于提高地质雷达信号分辨率的方法 |
| CN109581364A (zh) * | 2018-12-14 | 2019-04-05 | 江苏中路工程技术研究院有限公司 | 基于探地雷达的桥梁搭板脱空定量化识别与自动标识方法 |
| CN109685011A (zh) * | 2018-12-25 | 2019-04-26 | 北京华航无线电测量研究所 | 一种基于深度学习的地下管线检测识别方法 |
| CN111175761A (zh) * | 2019-11-19 | 2020-05-19 | 南京工程学院 | 水下机器人定位声纳数据的配准方法 |
| CN111323774B (zh) * | 2020-03-30 | 2022-06-14 | 华南农业大学 | 一种采用几何圆柱探测模型在探地雷达图谱中提取双曲线信号的方法 |
| CN112924956B (zh) * | 2021-02-05 | 2022-05-13 | 中国人民解放军国防科技大学 | 一种降低低频噪声级的背景均衡方法 |
| CN113688692A (zh) * | 2021-07-29 | 2021-11-23 | 浙江大学 | 基于探地雷达时频尺度特征的供水管道漏损检测方法 |
Citations (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| EP2037293A1 (en) * | 2006-12-21 | 2009-03-18 | Galileo Avionica S.p.A. | Multiple-target radar recognition method and apparatus |
| CN102737248A (zh) * | 2012-06-21 | 2012-10-17 | 河南工业大学 | 复杂路况下的车道线特征点提取方法与装置 |
| CN103902985A (zh) * | 2014-04-15 | 2014-07-02 | 安徽工程大学 | 一种基于roi的强鲁棒性实时车道侦测算法 |
| CN104020495A (zh) * | 2014-06-24 | 2014-09-03 | 中国矿业大学(北京) | 一种基于探地雷达的地下管线参数自识别方法 |
| CN105005042A (zh) * | 2015-07-27 | 2015-10-28 | 河南工业大学 | 一种探地雷达地下目标定位方法 |
-
2015
- 2015-12-16 CN CN201510942842.3A patent/CN105353373B/zh active Active
Patent Citations (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| EP2037293A1 (en) * | 2006-12-21 | 2009-03-18 | Galileo Avionica S.p.A. | Multiple-target radar recognition method and apparatus |
| CN102737248A (zh) * | 2012-06-21 | 2012-10-17 | 河南工业大学 | 复杂路况下的车道线特征点提取方法与装置 |
| CN103902985A (zh) * | 2014-04-15 | 2014-07-02 | 安徽工程大学 | 一种基于roi的强鲁棒性实时车道侦测算法 |
| CN104020495A (zh) * | 2014-06-24 | 2014-09-03 | 中国矿业大学(北京) | 一种基于探地雷达的地下管线参数自识别方法 |
| CN105005042A (zh) * | 2015-07-27 | 2015-10-28 | 河南工业大学 | 一种探地雷达地下目标定位方法 |
Non-Patent Citations (3)
| Title |
|---|
| 利用互相关和Hough变化快速检测探地雷达目标;王建等;《电子与信息学报》;20130531;第35卷(第5期);第1156-1162页 * |
| 基于广义Hough变换的探地雷达图像地下管线提取方法;孙伟等;《测绘科学技术学报》;20131231;第30卷(第3期);第251-259页 * |
| 浅地层探地雷达波速测量方法的研究;孔令讲;《电子学报》;20020930;第30卷(第9期);第1331页 * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN105353373A (zh) | 2016-02-24 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN105353373B (zh) | 一种基于Hough变换探地雷达目标提取方法和装置 | |
| CN105005042B (zh) | 一种探地雷达地下目标定位方法 | |
| CN109521422B (zh) | 一种基于雷达信号的多目标生命探测方法及探测雷达 | |
| CN104020495B (zh) | 一种基于探地雷达的地下管线参数自识别方法 | |
| Martorella | Novel approach for ISAR image cross-range scaling | |
| CN108387896B (zh) | 一种基于探地雷达回波数据的自动汇聚成像方法 | |
| US11900630B2 (en) | Method for detecting leakage of water supply pipe based on ground-penetrating radar three-dimensional image attribute analysis | |
| CN109343022B (zh) | 估测层间土壤含水量的方法 | |
| CN106443775B (zh) | 高分辨率转换波裂缝预测方法 | |
| CN106022339B (zh) | 一种复垦土地浅埋地埋管深度的提取方法 | |
| CN109085572B (zh) | 隧道内基于多径利用毫米波雷达的运动目标跟踪方法 | |
| CN110133643B (zh) | 植物根系探测方法及装置 | |
| CN106526678B (zh) | 一种反射声波测井的波场分离方法及装置 | |
| CN107329137B (zh) | 基于box样条小波框架的探地雷达地下目标定位方法及系统 | |
| CN109581317B (zh) | 一种基于回波峰值匹配的角落目标定位方法 | |
| CN106019254A (zh) | 一种uwb冲击生物雷达多人体目标距离向分离辨识方法 | |
| CN107220628B (zh) | 红外干扰源检测的方法 | |
| CN108802725A (zh) | 一种浅层穿透雷达合成孔径成像方法 | |
| CN111679275A (zh) | 基于探地雷达的地下管线识别方法 | |
| Qin et al. | Bidimensional empirical mode decomposition method for image processing in sensing system | |
| Luo et al. | Accurate tree roots positioning and sizing over undulated ground surfaces by common offset GPR measurements | |
| CN110261905A (zh) | 基于倾角控制的复值相干微断层识别方法 | |
| CN113780361A (zh) | 一种基于2.5d-cnn算法的三维探地雷达图像地下管线识别方法 | |
| CN113985487A (zh) | 一种基于三维探地雷达的地下埋藏物三维渲染方法及系统 | |
| CN107255815A (zh) | 一种基于双基地散射中心时频特征的目标曲面重构方法 |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| C06 | Publication | ||
| PB01 | Publication | ||
| C10 | Entry into substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |