Fast detection of GPR objects with cross correlation and Hough transform.
In most Ground Penetrating Radar (GPR) applications, object detection is a crucial but difficult problem. An efficient approach is typically based on detection of a certain spatial distribution for reflected energy . The hyperbola is a simple but effective distribution [2-6]. Based on a hyperbolic distribution, the Hough transform has successfully been used to detect GPR objects by detecting hyperbolas in a B-scan [2, 3]. However, the multiple parameters for a hyperbola yield a low calculation efficiency for Hough-based detection methods, which results in their limited practical application. Certain improved Hough transforms, such as the fast  and probabilistic Hough transforms , have reduced the parameter space using thresholds. However, the thresholds must be determined on a per problem basis. Furthermore, a binary image of edges is used as input, which is a difficult task for GPR images. The same question exists in the paper .
In this paper, we propose a fast Hough transform-based detection algorithm. The object reflection is modeled as a three-parameter hyperbola. To enhance computation efficiency, we first estimate two parameters for the reflection hyperbola; thus, only a 1D Hough transform is necessary to detect the hyperbola. As for parameter estimation, we utilized the correlation between adjacent 1D reflections (A-scans) to extract a hyperbola, which can be used to estimate two parameters for the reflection hyperbola, and the extracted curve is weighted using the reflected energy and fitted to derive the two parameters. The dimension reduction for the Hough transform will likely reduce the computational load.
2. REFLECTION MODEL OF A BURIED OBJECT
In a B-scan, the reflected energy of a small subsurface object is approximately distributed as a hyperbola, as shown in Figure 1. To simplify the problem, the GPR system is assumed monostatic, and the underground medium is assumed homogeneous. When the ground-coupled antenna scans an object linearly, the antenna position x and corresponding echo time delay t approximately satisfy a three-parameter hyperbola equation [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]:
[[t.sup.2]/[t.sup.2.sub.0]] - [4[(x - [x.sub.0]).sup.2]/[[upsilon].sup.2][t.sup.2.sub.0]] = 1 (1)
where [upsilon] is the velocity of wave propagation underground, [x.sub.0] the horizontal position of the object, and [t.sub.0] the time delay for the echo when the antenna is above the object. The apex ([x.sub.0], [t.sub.0]) [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] indicates the object position. Along the reflection hyperbola, the reflected energy is weaker as [absolute value of x - [x.sub.0]] increases.
3. DETECTION ALGORITHM
3.1. Extract Reflection Hyperbola by Cross Correlation
Correlation receiver is a widely used method of detecting objects in free space. In this method, the received signal y(t) is modeled as follows
y(t) = [1/a] s(t - [tau]) + w(t) (2)
where s(t) is the transmitted signal. [tau] and 1/a (a > 1) are the time delay for an object reflection and amplitude attenuation factor, respectively, and w(t) represents noise. The cross correlation between the transmitted and received signals is the following
r(t) = [summation over (u)] y(u)s(u - t) (3)
r([??]) > [gamma], [??] = arg max r(t) (4)
then an object is detected. [gamma] is a threshold derived from the given false alarm rate (FAR).
However, the correlation receiver does not operate well for a GPR because the waveform of the received GPR signal is distorted by the dispersive media underground, and the received signal model described by (2) is not valid. Nevertheless, for GPR, two adjacent A-scans are strongly correlated. Following the notion of a correlation receiver, we modeled the GPR received signal as follows
[A.sub.n](t) = [1/a] [A.sub.n-1](t - [DELTA][[tau].sub.n]) + w(t), n = 2, 3, ... (5)
where [A.sub.n](t) represents the nth A-scan, 1/a (a > 1) the amplitude attenuation factor, w(t) the noise, and [DELTA][[tau].sub.n] the relative time delay for [A.sub.n] (t) relative to [A.sub.n-1](t). Where [[tau].sub.n] denotes the object reflection time delay of [A.sub.n](t), [DELTA][[tau].sub.n] is defined as follows
[DELTA][[tau].sub.n] = [[tau].sub.n] - [[tau].sub.n-1], n = 2, 3, ... (6)
The cross correlation between [A.sub.n](t) and [A.sub.n-1](t) is as follows
[R.sub.n](t) = [summation over (u)] [A.sub.n](u) [A.sub.n-1](u - t) (7)
After the direct wave is suppressed by preprocessing , [DELTA][[tau].sub.n] can be estimated using the following equation
[DELTA][[tau].sub.n] = arg max [R.sub.n](t) (8)
According to previous studies (6), [[tau].sub.n] can be derived recursively as
[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (9)
The point pairs ([x.sub.n], [[tau].sub.n]) compose the reflection hyperbola [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII], [x.sub.n] is the antenna position where [A.sub.n](t) is collected. However, we cannot calculate [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] because it is difficult to estimate the initial value [[tau].sub.1]. Nevertheless, we can estimate two parameters for [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII], [x.sub.0] and [upsilon], from ([x.sub.n], [[tau].sub.n]) with any [t.sub.1], which is why ([x.sub.n], [[tau].sub.n]) was extracted.
For any given initial value [[tau]'.sub.1], the following applies
[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (10)
Let [DELTA][[tau].sub.1] = [[tau].sub.1] - [[tau]'.sub.1], then
[[tau]'.sub.n] = [[tau].sub.n] - [DELTA][[tau].sub.1], n = 1, 2, ... (11)
Therefore, the equation for the hyperbola ([x.sub.n], [[tau]'.sub.n]) is [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII], which indicates that [x.sub.0] and [upsilon] can be estimated by fitting ([x.sub.n], [[tau]'.sub.n]). For simplicity, [[tau].sub.1] is set to zero.
Figure 2 shows an example of reflection hyperbola extraction. Figure 2(a) is a simulated B-scan, and Figure 2(b) is the extraction result where [[tau].sub.1] = 0. In Figure 2(c), we placed the extracted hyperbola and reflection hyperbola with two apexes aligned along direction T. The two hyperbolas are the same except for a shift in the direction T, [x.sub.0] and [upsilon] can be estimated from the extracted hyperbola.
In practical applications, noise will deteriorate the estimation of [DELTA][[tau].sub.n], but [x.sub.0] and [upsilon] can still be correctly estimated. Like Figure 2, Figure 3 illustrates extraction in a noisy B-scan. Now, however, the extracted curve can be divided into three segments. The middle and smooth segments comprise the valid portion, which is slightly affected by noise. The side and rough segments are the invalid portion, which is seriously affected by noise. The comparison in Figure 3(c) shows that the valid portion of the extracted curve is consistent with the reflection hyperbola, which can be explained using two observations. First, the signal to noise ratio (SNR) of the middle portion in the reflection hyperbola is sufficiently high to estimate [DELTA][[tau].sub.n]. Second, the extraction is recursive, thus, the left invalid portion only shifts the portion extracted thereafter in the direction T, similar to the initial value [[tau].sub.1], and the right invalid portion has no impact on the pre-extracted portion. Therefore, if the SNR in the middle portion of a reflection hyperbola is sufficiently high, which is usually true, [x.sub.0] and [upsilon] can be estimated by fitting the valid portion of the extracted hyperbola.
3.2. Estimating Two Parameters for a Reflection Hyperbola
As a quadratic curve, the extracted hyperbola can be easily fitted using the least squares method (LS) to derive the two parameters we want. However, if the extracted curve is directly fitted, the invalid portions arisen by the noise will interfere in the parameter estimation, which is illustrated in Figure 4. To estimate [x.sub.0] and [upsilon], we only need to fit the valid portion but the entire curve. Thus, we use weighted fitting to suppress the invalid portion.
In fact, the cross correlation value [R.sub.n]([DELTA][[tau].sub.n]) can be used as the weight because [R.sub.n]([DELTA][[tau].sub.n] is weaker when [absolute value of [x.sub.n] - [x.sub.0]] increases, and [R.sub.n]([DELTA][[tau].sub.n]) in the invalid portion is much smaller than the valid portion. For example, Figure 5 shows the cross correlation value for the extracted hyperbola depicted in Figure 3(b), and the weighted fitting result is in Figure 4. We find that the weighted fitting suppresses the invalid portion of the extracted hyperbola and is consistent with the valid portion. Further, [x.sub.0] and [upsilon] can be estimated by fitting the extracted hyperbola using [R.sub.n]([DELTA][[tau].sub.n]) as the weight.
3.3. Detecting and Localizing Objects Using a 1D Hough Transform
For the traditional Hough-based GPR object detection method, the 3D Hough transform is necessary, as follows:
[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (12)
where B(x,t) is the B-scan data. If a hyperbola exists in a B -scan, a peak is generated in [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]. Further, objects are detected through detecting peaks in [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII], while object localization can use the following equation
[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (13)
In traditional Hough-based methods, a Hough transform, peak detection and localization are all considered in a 3D parameter space [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII], which consumes a long computational time.
Because we estimate two parameters for the reflection hyperbola, unlike the traditional detection method, only a 1D Hough transform is necessary, as follows
[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (14)
Similarly, the peak in H([[??].sub.0]) indicates the existence of the hyperbola [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]. Further, object localization can be implemented by deriving [t.sub.0], the rest parameter of the hyperbola model, from the following equation
[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (15)
Compared with (12) and (13), dimensionality reduction in (14) and (15) distinctly reduces the computational cost of our method.
4. EXPERIMENTAL RESULTS AND ANALYSIS
4.1. Experiment Data
The experimental data were collected by Radar Eye, which is a pulse GPR system developed by our group. The Radar Eye antennas are a pair of TEM rectangular horns, and the transmitting signal is a bipolar pulse with the central frequency 1 GHz. Figure 6 shows the experimental setup. Small objects are buried in dry sand, and the antennas are 5 cm above the sand. We collected 150 B-scans, including 80 metallic small objects, 20 plastic objects, and 50 scans of random clutter (responses produced by random noise or an unknown mechanism). The size of each object was less than 10cm x 10cm x 10 cm, and the burial depth was between 10 cm to 20 cm. Figure 7 shows 3 B-scan samples, which are collected from a metallic object, a plastic object, and clutter.
There are many methods for detecting small underground objects. Methods such as Hidden Markov models (HMM) , spectral confidence feature (SCF), edge histogram discrimination (EHD) , and polynomial fitting (PF)  have been successful at landmine detection. Herein, we compare the proposed method to HMM, PF and the method based on standard Hough transform (SHT). The programs were coded and executed in Matlab 2009a using a computer with a 3.06 GHz Pentium[R] 4 CPU.
Figure 8 shows the ROC results. The probability of detection (PD) is defined as the number of targets detected divided by the total number of targets (i.e., 100). The FAR is defined as the ratio for the number of detected clutters to the total clutters (i.e., 50). The results from HMM, SHT and our method can produce 90% PD at 8% FAR. In this experiment, the PF performs comparatively worse, but 90% PD can be generated at 10% FAR. When the FAR is less than 6%, the HMM performs better than the other algorithms. Overall, the detection ability of our method is encouraging.
The average times for a B-scan are compared in Table 1. The B-scans are the same size, with 2750 time samples and 100 position samples. It shows that the computational time for our method and PF is approximately 2% of the HMM and SHT time, which indicates that our method is highly efficient.
Our method's advantageous calculation efficiency can be proven theoretically. The size of the B-scan is assumed M xN, where M and N are the sample numbers for time and space, respectively. The dominant calculation for the traditional Hough-based method uses a 3D Hough transform. The times used by multiplication and summation processes in a 3D Hough transform produce are both O(M x N2 x Nv), where Nv is the sample number of [upsilon]. In our method, the calculation primarily comes from 3 steps, computing cross correlation, fitting hyperbola and a 1D Hough transform. The times used by the multiplication and summation processes during the cross correlation computation are both O(M x N x L), where L is the correlation length. L is typically short because the relative time delay for two adjacent received signals is small. Thus, O(M x N x L) is close to O(M x N). To fit a quadratic curve with 3 parameters, both multiplication and summation use O(N) time. The time for multiplication and summation in a 1D Hough transform are O(M x N). Table 2 compares the computational cost of two methods. The time used by a traditional method is proportional to M x [N.sup.2] x [N.sub.[upsilon]], while our method is M x N. Thus, our method is much faster than the traditional method.
In this paper, a fast method based on cross correlation and a hough transform is presented for GPR object detection. The method is compared with three other detection methods using filed data. The results show an encouraging detection ability and high computational efficiency. The high efficiency was also analyzed theoretically. The proposed algorithm is valuable for detecting small underground objects when high computational efficiency is required. However, in practical applications, the inhomogeneity of the medium underground and the unstable antenna height will degrade the detection performance.
[1.] Jol, H. M., Ground Penetrating Radar Theory and Applications, Elsevier Science, Amsterdam, 2009.
[2.] Carlotto, M. J., "Detecting buried mines in ground penetrating radar using a hough transform approach," Battlespace Digitization and Network-Centric Warfare II, 251-261, Orlando, 2002.
[3.] Chen, D., C. Huang, and Y. Su, "An integrated method of statistical method and hough transform for GPR target S detection and location," Acta Electronic Sinica, Vol. 32, No. 9, 1468-1471, 2004.
[4.] Simi, A., S. Bracciali, and G. Manacorda, "Hough transform based automatic pipe detection for array GPR: Algorithm development and on-site tests," IEEE Radar Conference, 1-6, Rome, 2008.
[5.] Birkenfeld, S., "Automatic detection of reflexion hyperbolas in GPR data with neural networks," World Automation Congress, 1-6, Kobe, 2010.
[6.] Liu, Y., M. Wang, and Q. Cai, "The target detection for GPR images based on curve fitting," 3rd International Congress on Image and Signal Processing, 2876-2879, Yantai, 2010.
[7.] Guil, N., J. Villalba, and E. L. Zapata, "A fast hough transform for segment detection," IEEE Trans. on Image Process., Vol. 4, No. 11, 1541-1548, 1995.
[8.] Kiryati, N., Y. Eldar, and A. M. Bruckstein, "A probabilistic hough transform," Pattern Recognition, Vol. 24, No. 4, 303-316, 1991.
[9.] Long, K., P. Liatsis, and N. Davidson, "Image processing of ground penetrating radar data for landmine detection," Proc. of SPIE, Vol. 6217, 62172R1-62172R12, 2006.
[10.] Hayashi, N. and M. Sato, "F-K filter designs to suppress direct waves for bistatic ground penetrating radar," IEEE Trans. on Geosci. Remote Sensing, Vol. 48, No. 3, 1433-1444, 2010.
[11.] Gader, P. D., M. Mystkowski, and Y. Zhao, "Landmine detection with ground penetrating radar using hidden Markov models," IEEE Trans. on Geosci. Remote Sensing, Vol. 39, 1231-1244, 2001.
[12.] Wilson, J. N., P. Gader, W.-H. Lee, H. Frigui, and K. C. Ho, "A large-scale systematic evaluation of algorithms using ground-penetrating radar for landmine detection and discrimination," IEEE Trans. on Geosci. Remote Sensing, Vol. 45, No. 8, 2560-2572, 2007.
[13.] Zhu, Q. and L. M. Collins, "Application of feature extraction methods for landmine detection using the Wichmann/Niitek ground-penetrating radar," IEEE Trans. on Geosci. Remote Sensing, Vol. 43, No. 1, 81-85, 2005.
Jian Wang * and Yi Su
School of Electronic Science and Engineering, National University of Defense Technology, Changsha 410073, China
* Corresponding author: Jian Wang (email@example.com).
Received 25 February 2013, Accepted 1 April 2013, Scheduled 2 April 2013
Table 1. Average time for a B-scan. Detection method Time (s) HMM 20.17 PF 0.31 SHT 23.34 Our method 0.35 Table 2. Computational cost comparison. Detection method Computational time Multiplication Traditional Hough-based method O(M x [N.sup.2] x [N.sub.[upsilon]]) Our method Cross correlation O(M x N) Hyperbola fitting O(N) 1D Hough transform O(M x N) Total O(M x N) + O(N) Detection method Computational time Summation Traditional Hough-based method O(M x [N.sup.2] x [N.sub.[upsilon]]) Our method Cross correlation O(M x N) Hyperbola fitting O(N) 1D Hough transform O(M x N) Total O(M x N) + O(N)
|Printer friendly Cite/link Email Feedback|
|Title Annotation:||ground penetrating radar|
|Author:||Wang, Jian; Su, Yi|
|Publication:||Progress In Electromagnetics Research C|
|Date:||May 1, 2013|
|Previous Article:||Cross-polarization reduction of E-shaped microstrip array using spiral-ring resonator.|
|Next Article:||Improved method of node and threshold selection in wavelet packet transform for UWB impulse radio signal denoising.|