## 2 Materials and Methods
### Multi-Frequency GPR Data Fusion  
#### 2.2.1 Data Fusion Method  
As shown in Figures 2 and 3, data fusion refers to integrating the multiple-frequency GPR signals into one representation to enhance the detection information and obtain more reliable and detailed results. Suppose there are \(m\) different frequencies of GPR signals distributed in the length, \(L\), sliding window and that the sliding window will move \(\Delta L\) along the time series axis at each step. In the sliding process, the sliding window will contain new GPR signals and remove the oldest GPR signal from the window. The overall optimization model of all GPR signals, \(X(t)\), within the sliding window can be calculated by minimizing the cost function of the Markov norm as  
\[X^{*}(t)=\underset{K}{\text{argmin}}\big{\{}\sum\lVert X(t,K)\rVert-\sum \lVert X(t,K-\Delta K)\rVert+\sum\lVert X(t,K+\Delta K)\rVert\big{\}} \tag{3}\]  
where \(X^{*}(t)\) represents the overall optimization model of all radar signals within the sliding window, \(K\), at time \(t\); \(X(t,K)\) denotes the GPR signal data model within the sliding window,\(L\), at time \(t\); \(X(t,K-\Delta K)\) and \(X(t,K+\Delta K)\) represent the GPR signal data model within the sliding windows, \(K-\Delta K\) and \(K+\Delta K\), at time \(t\), respectively; and \(K\) and \(\Delta K\) are the window length and the moving step, respectively. For convenience, assume that the GPR signals within sliding window, \(L\), can be divided into three categories: high-, moderate- and low-frequency GPR signals. Then, the weighted fusion function of each segment frequency signal in the sliding window can be expressed as  
\[\widetilde{X}(t)=\omega_{H}X^{*}_{H}(t)+\omega_{M}X^{*}_{M}(t)+\omega_{L}X^{*} _{L}(t) \tag{4}\]  
where \(\omega_{H}\), \(\omega_{M}\) and \(\omega_{L}\) denote the weighting factors of high-frequency, moderate-frequency and low-frequency wavelet signals in each segment of the sliding window, the values of \(\omega_{H}\), \(\omega_{M}\) and \(\omega_{L}\) are within 0 and 1, also, \(\omega_{H}+\omega_{M}+\omega_{L}=1\). Next, we use the wavelet transform method to calculate the time-varying weight factors of the GPR signals in every frequency band within the sliding window. Suppose that \((x,y)\) is the center point of the sliding window, \(L\); \(\sigma_{H}\), \(\sigma_{M}\) and \(\sigma_{L}\) denote the root mean square errors of high-, moderate- and low-frequency GPR signals within the sliding window, respectively, with variances of \(D(X_{H})=\sigma_{H}^{2}\), \(D(X_{M})=\sigma_{M}^{2}\) and \(D(X_{L})=\sigma_{L}^{2}\). To  
Figure 3: The weighting joint of the sliding window and wavelet transform. ensure the \(D(X)\) optimization, we should calculate the minimum value of \(D(X)\) under the constraint condition \(\omega_{H}+\omega_{M}+\omega_{L}=1\), which is expressed as  
\[\left\{ \begin{array}{l}\omega_{H}=\frac{\frac{\sum\limits_{i}^{T}\frac{y}{L} (\sigma_{H}^{2}(t)w_{Li}^{l}(i,j))}{\sigma_{H}^{2}(t)+\sigma_{M}^{2}(t)+\sigma _{L}^{2}(t)}}{\sigma_{H}^{2}(t)+\sigma_{M}^{2}(t)+\sigma_{L}^{2}(t)}\\ \omega_{M}=\frac{\sum\limits_{i}^{T}\frac{y}{L}(\sigma_{L}^{2}(t)w_{Li}^{l}(i,j))}{\sigma_{H}^{2}(t)+\sigma_{M}^{2}(t)+\sigma_{L}^{2}(t)}\\ \omega_{L}=\frac{\sum\limits_{i}^{T}\frac{y}{L}(\sigma_{L}^{2}(t)w_{Li}^{l}(i,j))}{\sigma_{H}^{2}(t)+\sigma_{M}^{2}(t)+\sigma_{L}^{2}(t)}\end{array} \right. \tag{5}\]  
where \(W_{H}^{l}(i,j)\), \(W_{M}^{l}(i,j)\) and \(W_{L}^{l}(i,j)\) denote the wavelet function of the high-, moderate- and low-frequency GPR signals within the wavelet decomposition scale, \(l\), and \(l\) is the wavelet decomposition scale. According to the above discussion, we can substitute Equations (3) and (5) into Equation (4) to obtain the multi-frequency GPR fusing result. The fusing step is summarized in Figure 3; it includes the following steps:  
1. Obtaining the preprocessing and spatial alignment for the input multi-frequency GPR images. 2.