The Joint UAV-Borne Magnetic Detection System and Cart-Mounted Time Domain Electromagnetic System for UXO Detection

Abstract

For unexploded O=ordnance (UXO) detection, individual technology cannot achieve the best detection performance. The new detection mode of joint magnetic and electromagnetic method has attracted more and more attention. In this paper, a newly developed joint detection system is introduced, a multi-rotor UAV-based magnetic system (UAVMAG) and a cart-based time-domain electromagnetic detection system (TDEM-Cart) are combined, and the cooperative processing of magnetic field and electromagnetic data is proposed. The result of the joint inversion fuses the feature vector retrieved from the magnetic field data and the feature vector inverted from the electromagnetic data, providing more accurate positioning results and richer information, which is favorable to locate and distinguish the UXO. Two field experiments are conducted, and the results show that when the joint system works in the full-coverage survey mode, both ferromagnetic and non-ferromagnetic metal targets can be detected, avoiding missed detections. In addition, when the joint system works in the cued survey mode, the detection efficiency is improved, the positioning accuracy of joint interpretation is less than 10 cm, and it shows satisfactory performance in the recognition of targets.

1. Introduction

Unexploded ordnance (UXO) is a product left over from a war, usually buried within 30 m underground [1,2]. With the expansion of human activities, potential unexploded ordnances pose a serious threat to the lives of residents, the living environment, and the economic development of society. Since the middle of the 19th century, countries around the world have been actively carrying out the detection and removal of UXO [3]. At present, multiple geophysical techniques including magnetic method [4,5,6], electromagnetic method [7,8] and ground penetrating radar (GPR) [9,10] are applied in the detection of UXO, and there are already a variety of prototypes and commercial products designed for the locating and identifying of UXO.

Magnetic survey is a kind of passive detection technique. Affected by the magnetization effect of the geomagnetic field, the ferromagnetic target produces an induced magnetic field, causing a local magnetic anomaly field [3]. Various types of scalar or vector magnetic field sensors are fixed on different carrier platforms, forming the handheld [11,12], helicopter-borne [13] and UAV-mounted magnetic measurement systems [14,15]. The drone-based magnetic detection system shows significant advantages in detection efficiency and spatial resolution of the magnetic data [16], which is a hot spot of current research.

Electromagnetic method is divided into two categories: Time Domain Electromagnetic Method (TDEM) and Frequency Domain Electromagnetic Method (FDEM), both of which consist of transmitting and receiving. The transmitter emits the primary field by controlling the current of the transmitting coil, and the receiver picks up the secondary field, which is an induced magnetic field produced by the metal target. Typical commercial products, such as EM61, EM61-MK2A [7], EM63-3D [8], and MPV-II [17], have been widely used in UXO detection, underground cable detection, archaeology, etc. Similarly, the UAV-based EM detection system [18,19,20] has attracted extensive interests, and the research and integration of the UAV-TEM system is more difficult than the UAV-magnetic system [21], but it is bound to play an important role in geological work with its advantages in low-altitude detection.

The detection and identification of UXO remains an ongoing challenge due to the limitations of terrain, the material of UXO and various clutters, a single detection technique cannot achieve the best detection performance. Magnetic detection has remarkable advantages in fast monitoring and large-scale survey, but it can only detect targets with high magnetic permeability. Moreover, it is difficult to quantitatively retrieve the shape and attitude information of the target based on the magnetic field data. In contrast, the electromagnetic method is able to detect metal targets including magnetic and non-magnetic, and the conductivity, magnetic permeability, shape, position, and attitude of the target can be inverted based on the rate of decay and spatial behavior of the secondary field [22,23,24]. However, the EM method is more complex and time consuming, and its detection distance is shallower than that of magnetic detection [23]. The new detection mode of combined magnetic and electromagnetic method has emerged as an interesting candidate [23,24,25,26]. Joint detection and inversion take advantage of the strengths of both techniques to help improve the detection efficiency and target recognition capabilities.

In terms of joint detection systems, a variety of combined modes based on vehicle, airborne or handheld have been developed and tested. GEO-CENTERS, Inc., in Masschusetts, U.S.A and the U. S. Army Corps of Engineers developed a vehicle-mounted joint detection system, one Geonics EM61 pulse induction sensor and multiple Geometrics 822A magnetometers were simultaneously towed by a non-magnetic cart [24]. A sophisticated sampling scheme was proposed to eliminate the electromagnetic interference of EM61 to the magnetometer, the magnetic data and EM61 data are synchronized in space based on the position information measured by GPS. The UXO detection tests at the Calibration Grid and the Blind Test Grid have verified that the detection efficiency of the combined vehicle-towed system has been improved significantly and the joint detection system performed well when the two targets are closely distributed. Jacob Sheehan et al. [25] designed a joint detection system based on the helicopter flying platform for the large-scale detection of UXO. The VG-22 and TEM8 systems were, respectively, integrated on the Bell 206 Long Ranger helicopter. The VG-22 system collected the vertical gradient magnetic field, seven vertical gradient pods are integrated on the forearm of the aircraft, and each gradient pod contains two cesium vapor magnetometers vertically separated by 50 cm. The TEM8 system collected the secondary field through eight receiving coils fixed under the aircraft. Experimental results indicated that the combined detection could reduce the false alarm rate. Kenji Okazaki et al. used an airborne electromagnetic system (GREATEM) and a helicopter-borne magnetic system (HMS) to delineate the geological structures of long tunnel construction sites planned in the Otoineppu area of Hokkaido, Japan [26]. The GREATEM survey and HMS survey were carried out successively, obtaining the resistivity distribution map and magnetic anomaly map of the same test area. This combined survey provides richer information for the effective interpretation of underground structures.

In terms of cooperative interpretation of data, different methods for joint processing of magnetic field and electromagnetic data have been proposed. Jacob Sheehan [24] used the normalization and DC removal method to ensure that the aerial magnetic data and the aerial electromagnetic data are comparable, and the detection results of the two methods are merged together, achieving the “Decision Fusion” [27]. Leonard R. Pasion et al. proposed a cooperative inversion method, where information from the inversion of magnetic data is used as a constraint for inverting TDEM data [23]. The cooperative inversion consists of three steps: the first is the inversion of the magnetic data, the second is the inversion of TEM data with location constraints, and finally to estimate the shape and size of target from the magnetic data. The field tests at Yuma Proving Ground have confirmed that this cooperative interpretation performed better than the single inversion method in target classification and recognition.

Due to the limitations of mutual interference between systems, data fusion methods and the cost, the existing joint detection system and interpretation method are still in the early stage and not yet mature. For the helicopter-based joint detection systems, the spatial resolution of aerial data is much greater than that of terrestrial data, and a fine survey in small areas cannot be achieved. At the same time, low-altitude flight by a helicopter places high demand on pilots, which is a severe challenge. For the vehicle-mounted joint detection system, it is difficult to completely eliminate the influence of the primary field emitted by the electromagnetic detection system on the magnetometer, and the detection efficiency of the vehicle-mounted system is relatively low. In contrast, UAV flight platform is a good compromise to improve detection efficiency and spatial resolution and is able to fly autonomously at a low altitude.

In this paper, a newly developed joint detection system is introduced, a multi-rotor UAV-based magnetic detection system (UAVMAG) and a cart-based time-domain electromagnetic detection system (TDEM-Cart) are combined; and the complete joint processing flow is proposed for the detection and discrimination of UXO. Two detection modes are designed, one is full-coverage survey mode, and the other is cued survey mode. In the full-coverage survey mode, the same test area is scanned by UAVMAG and TDEM-Cart, respectively, ensuring that non-magnetic metal targets will not be missed. In the cued survey mode, UAVMAG system first quickly scans the entire test area, and a series of data processing methods are used to delineate the abnormal areas where the target might be, providing “clues” for the subsequent detection by TDEM-Cart. Next, the local abnormal area is finely surveyed through the TDEM-Cart, instead of covering the entire test area, thereby improving the detection efficiency. Furthermore, the cooperative data interpretation method is presented, magnetic data and electromagnetic data (EM data) are spatially synchronized by the Real-Time Kinematic Global Positioning System (RTK GPS) system, and the synchronization deviation is cm level. Before the quantitative interpretation of data, some data processing methods are adopted to improve the quality of magnetic field and EM data. The global optimization algorithm differential evolution (DE) is employed in the inversion of magnetic data, the local optimization method Levenberg–Marquardt (LM) is used in the inversion of EM data, and the feature vector retrieved from the magnetic data is set as the initial value of the electromagnetic data inversion. The joint interpretation can avoid falling into the local optimal solution, the positioning accuracy of the target is higher than the magnetic inversion result, and the information of the merged feature vector is richer, which is helpful to distinguish and recognize the target.

The paper is organized as follows. Section 2 introduces the composition of the joint detection system, including the basic structure of UAVMAG and TDEM-Cart. Section 3 proposes the cooperative interpretation of magnetic data and EM data for target positioning and recognition. Two experiments are conducted, and the results are discussed in Section 4. Section 5 gives the discussion and remarks of the joint detection system and processing method. Section 6 provides the conclusions of this research.

2. System

The joint detection system combines the UAV-borne magnetic measurement system (UAVMAG) and the cart-mounted time domain electromagnetic detection system (TDEM-Cart). In this paper, UAVMAG and TDEM-Cart do not work at the same time, so there is no need to consider the interference between the two systems. Both systems are equipped with the RTK GPS, and a unified spatial coordinate system is established to ensure the spatial synchronization of magnetic field data and electromagnetic data.

2.1. UAVMAG System

UAVMAG is a passive detection system, and magnetic field sensor is the main component of the system. Figure 1 shows the structure of UAVMAG system. Two optically pumped magnetometers are fixed under the drone through a vertical rigid rod, and the connector between the rod and the drone can be rotated within 90 degrees, ensuring that the rigid rod is switched from horizontal to vertical, achieving the automatic take-off and landing for UAVMAG.

The flight platform is a commercial eight-rotor drone MG1 produced by SZ DJI Technology Co., Ltd., China. The optically pumped magnetometer is a scalar magnetic field sensor CAS-L3 developed by The Chinese Academy of Sciences. The noise sensitivity of CAS-L3 is 0.6 pTrms√Hz⁄@1 Hz in a shielded room, absolute accuracy is less than 2.5 nT (nano Tesla). UAVMAG is also equipped with RTK GPS T300, which is manufactured by the Shanghai Sinan Satellite Navigation Technology Co., Ltd., providing positioning and timing service, the positioning accuracy is cm level, the accuracy of pps (pulse per second) [28] is 20 ns.

Two magnetometers and RTK GPS are synchronized by the pps signal. The acquisition system collects and stores the magnetic field and position information, and the sampling rate is 160 Hz. As a result, the magnetic survey throughout the test area could be realized using this UAVMAG system with a maximum flying speed up to 12 m/s. One line per 7.5 cm could be collected when the drone reaches its maximum flying speed.

In this paper, the UAVMAG system is a part of the joint detection system, which is used to collect spatial magnetic field data at different elevations, achieving the fast detection of large-scale areas.

2.2. TDEM-Cart System

TDEM-Cart is an active detection system, including two parts: transmitting system and receiving system. Figure 2 shows TDEM-Cart system and its electronic block diagram. The TDEM-Cart system consists of a transmitting coil and three receiving coils, RTK GPS, and the host module. All coils are mounted on a non-magnetic cart, three receiving coils are fixed along the vertical direction, with an interval of 27.5 cm. The host module contains H-bridge circuit, low noise amplifier (LNA), band pass filter (BPF), analog-to-digital converter (ADC), field programmable gate array (FPGA), and synchronous dynamic random-access memory (SDRAM), realizing the waveform design of the transmitting system and the signal acquisition and synchronization of the receiving system.

The transmitting coil is a 1 m × 1 m square coil, an H-bridge circuit is used to generate a bipolar pulse current, and FPGA is employed to achieve data transmission, timing control, frequency adjustment and the off time of the transmission current. Finally, a bipolar square wave with a 50% duty cycle is generated, the peak current of the emission current is 8A, the base frequency is adjustable from 3.125 Hz to 50 Hz, a simulated emission current waveform diagram is also given in Figure 2.

The dimensions of the three receiving coils (labeled as No.1, No.2 and No.3 receiving coil) are the same, all of which are 0.5 m × 0.5 m square coils, and the receiving system includes LNA, BPF, ADC, FPGA and SDRAM modules, as illustrated in Figure 2. Three-channel induced voltages are digitized at a sample rate of 156 kHz with 24-bit resolution, and RTK GPS T300 provides location and time data. The timing control, data synchronization and transmission are completed by FPGA, and thus the electromagnetic data, position and time data are synchronized and stored in the SDRAM. The schematic diagram of simulated single channel received voltage signal is also demonstrated in Figure 2.

In this paper, the TDEM-Cart system is a part of the joint detection system. The transmitting system of TDEM-Cart transmits a bipolar square wave, generating the primary field, and the receiving system of TDEM-Cart receives the secondary field generated by the target during the off time. For the location and discrimination of UXO, only the EM data of No.1 receiving coil are employed in subsequent processing.

3. Methods

The joint data processing flow is divided into four parts, the first part is the magnetic field signal processing, the second part is the inversion of magnetic data, the third part is the EM data processing, and the fourth is the inversion of EM data. The signal processing algorithm is executed before the inversion operation, which is used to extract the magnetic anomaly field or the secondary field signal related to the target, eliminating the noise signal and improving the data quality and signal-to-noise ratio (SNR). The fused feature vector containing the location, shape and posture information of the target is retrieved from the inversion results to locate and distinguish the target. The four parts are explained in detail below.

3.1. Magnetic Field Signal Processing

The measured result of the magnetic field sensor fixed on the UAV is the superposition of the geomagnetic field, the interference magnetic field generated by the flying platform, environmental noise, and the magnetic anomaly signal related to the target. What we are interested in is the magnetic anomaly signal; therefore, it is necessary to separate the magnetic anomaly signal from the measured magnetic field signal through a variety of signal processing methods.

3.1.1. Lowpass Filtering

In frequency domain, the energy of magnetic anomaly signals generated by a static metal target is mainly concentrated in low frequencies, lower than 1 Hz [29], and low-pass filter is able to remove the high-frequency noise, such as the AC magnetic field generated by the power grid. A digital low-pass filter with the cut-off frequency of 3 Hz is designed, cancelling the high-frequency noise from the measured magnetic field signal. Figure 3 displays the frequency spectrum of the magnetic field signal before and after filtering.

3.1.2. DC Removal

The filtered magnetic field signal consists of the geomagnetic field, UAV interference magnetic field and the magnetic anomaly signal. In a small area, the geomagnetic field is approximately stable and spatially uniform in a short period of time [30]. Here, the geomagnetic field is eliminated by removing the DC component.

The coordinate system o-xyz is established; x points to the east, y points to the north, and the z-axis is vertically upward, while the flying start point of the drone is used as the origin of the coordinate. The UAVMAG system flies a line in the north–south direction, flying speed of the drone is 2 m/s, and the sampling rate is 160 Hz. The raw magnetic field data is processed by a lowpass filter in Section 3.1.1, as shown in Figure 4a. Figure 4 compares the magnetic field before and after removing the DC component, and the geomagnetic field is eliminated in Figure 4b.

3.1.3. The Elimination of UAV Interference Field

The drone is not non-magnetic, partial components or structures are made of ferrous material, which will produce an interference magnetic field. The platform interference magnetic field is superimposed on the magnetic anomaly signal, deteriorating the quality of the magnetic data. The interference signals measured by two closely spaced sensors are correlated, and the interference signal is not correlated with the magnetic anomaly field; therefore, the magnetic anomaly field can be extracted from the total magnetic field based on the signal correlation model, and the details are described in Reference [14], which is provided in Supplementary Materials.

Figure 5a,b demonstrate the measured results of magnetometer 1 and magnetometer 2, low-pass filtering and DC removal have been executed, and the extracted magnetic anomaly field after removing the UAV interference field is shown in Figure 5c.

Magnetometer 1 is closer to the drone than magnetometer 2, and the interference field measured by magnetometer 1 is larger, the root mean square (rms) errors of the two measured magnetic signals are 5.485 nT and 2.663 nT. The magnetic anomaly signal is separated from the measured results according to the method in Supplementary Materials, the platform interference magnetic field is significantly reduced, rms is reduced to 1.077 nT.

3.2. The Inversion of Magnetic Field

For the ferromagnetic target in the shallow subsurface, when the distance between the sensor and the target is more than three times the maximum size of the target, the induced magnetic field generated by the target can be approximated as a point-shaped magnetic dipole. The dipole model is a point model, ignoring the shape and posture information of the target, and it can be characterized by six parameters including three-dimensional position vector and three-dimensional magnetic moment vector. The position vector provides the information about the location of the target, and the magnetic moment vector is related to the posture, shape, and material of the target.

A magnetic dipole is located at $r_0=\left(r_{0x},r_{0y},r_{0z}\right)$, the magnetic moment is $m=\left(m_x, m_y, m_z\right)$, and the magnetic field $B_{dipole}$ at the observation point $r_p$ is:

$$B_{dipole}\left(r,m\right)=\frac{{\mu }_0}{4\pi }\left(\frac{3\left(m\cdot r\right)r}{r^5}-\frac{m}{r^3}\right)$$

(1)

${\mu }_0$ is the permeability in the vacuum, and equals to $4\pi \times {10}^{-7}$ H/m, $r$ is the distance vector from the observation point $r_p$ to the dipole $r_0$, and $r$ is its modulus.

Assuming that there is no measurement noise, the output of the optically pumped magnetometer $B_{total}$ is the superposition of the geomagnetic field and the projection of the induced magnetic field in the direction of the geomagnetic field [31]:

$$B_{total}=\left|B_{earth}\right|+B_{dipole}·\frac{B_{earth}}{\left|B_{earth}\right|}$$

(2)

$B_{earth}$ is the geomagnetic field, $\left|B_{earth}\right|$ is its magnitude.

Consequently, the magnetic anomaly field $\Delta T$ is related to the magnetic moment and position of the dipole, written as:

$$\Delta T=f\left(r,m\right)=B_{dipole}·\frac{B_{earth}}{\left|B_{earth}\right|}$$

(3)

The inversion problem is to invert the position $r_0$ and magnetic moment $m$ based on the magnetic anomaly signal ($\Delta T_i$, i = 1, 2, …, N) at multiple observation points ($r_p$, i = 1, 2, …, N), six unknown parameters form the feature vector $x=\left(r_{0x},r_{0y},r_{0z},m_x,m_y,m_z\right)$. The differential evolution (DE) algorithm is adopted to estimate x, which is a classic heuristic random search algorithm [32,33]. DE is a global optimization strategy, suitable for solving noisy, multi-objective optimization problem.

The objective function is established by minimizing the data misfit between predicted magnetic anomaly field and observed magnetic anomaly field $\Delta T_i$. The predicted magnetic anomaly field $f_i\left(x\right)$ of the ferromagnetic target is based on Equation (3), so the objection function is expressed as follows:

$$\arg \underset{x}{\min }{\displaystyle \sum _{i=1}^N{\left(f_i\left(x\right)-\Delta T_i\right)}^2}$$

(4)

The inversion flow of magnetic data is listed in Algorithm 1. The feature vector x inverted from the magnetic data provides the position and magnetic moment of the buried target.

Algorithm 1. The Inversion of Magnetic Field.

Input: Magnetic anomaly field $\Delta T_i$; The range of feature vector: $x_{\min }^j$, $x_{\max }^j$ (j = 1,2,..,6) Population: M; Dimension: D; Generation: T Output: The location and magnetic moment vector of the target: $\widehat{x}$ 1  t$\leftarrow$1 (Initialization); 2  for i = 1 to M do 3 for j = 1 to D do 4  $x_{i,t}^j=x_{\min }^j+rand\left(0,1\right)\ast \left(x_{\max }^j-x_{\min }^j\right)$; 5 end 6  end 7  while $\left(t\le T\right)$ or $\left|f\left(\widehat{x}\right)\right|\ge \epsilon$ do 8 for i = 1 to M do 9  for j = 1 to D do 10   $v_{i,t}^j=Mutation\left(x_{i,t}^j\right)$; 11   $u_{i,t}^j=Crossover\left(x_{i,t}^j,v_{i,t}^j\right)$; 12  end 13  if $f\left(u_{i,t}\right)<f\left(x_{i,t}\right)$ then 14   $x_{i,t}\leftarrow u_{i,t}$; 15   if $f\left(x_{i,t}\right)<f\left(\widehat{x}\right)$ then 16    $\widehat{x}\leftarrow x_{i,t}$; 17   end 18  else 19   $x_{i,t}\leftarrow x_{i,t}$; 20  end 21 end 22 t$\leftarrow$ t + 1; 23  end 24  return feature vector $\widehat{x}$ including location and magnetic moment

3.3. EM Data Processing

TDEM-Cart system receives the secondary field generated by the target and converts it into the voltage signal through the induction coil. Its working principle is Faraday’s law of induction, and the working frequency band is very wide. The electromagnetic noise in the ambient environment is overlapped with the response of the target, the SNR of the secondary field will be reduced, degrading the detection performance. Therefore, some signal processing methods are utilized to eliminate the environmental noise and improve the SNR.

3.3.1. Signal Stacking

The emission current of the TDEM-Cart is a bipolar square wave, so the corresponding secondary field is also alternately positive and negative. Assuming that the received data has n ($n\ge 1$) periods, there are a total of 2n positive and negative response signals. The negative polarity response is reversed and superimposed with the positive polarity response:

$$V_s\left(t\right)=V_i\left(t\right)\ast s_n\left(t\right)$$

(5)

$V_i\left(t\right)$ and $V_s\left(t\right)$ are the signals before and after stacking, and $s_n\left(t\right)$ is the stacking function, expressed as:

$$s_n\left(t\right)=\frac{1}{2n}\left\{\delta \left(t\right)-\delta \left(t-\frac{T}{2}\right)+\delta \left(t-T\right)-\cdots +{\left(-1\right)}^{2n-1}\delta \left[t-\left(2n-1\right)\frac{T}{2}\right]\right\}=\frac{1}{2n}{\displaystyle \sum _{m=0}^{2n-1}{\left(-1\right)}^m\delta \left(t-m\frac{T}{2}\right)}$$

(6)

n denotes the number of the stacking period, T is the period of the emission current, $\delta \left(t\right)$ is the Dirichlet function.

To understand the Equations (5) and (6), the Fourier functions of the two stacking functions with n = 1 and n = 8 are demonstrated in Figure 6a. It shows that the stacking procedure eventually renders a narrow-band, match-filter-type band pass filter in time domain, and the larger the number of periods n, the narrower the bandwidth of the filter. Figure 6b demonstrates the same result in time domain, three response signals with different stacking cycle n = 0, 1, 8 are compared. As the number of superimposed cycles increases, the noise suppression is enhanced, especially in the late stage of the signal.

3.3.2. Windowing in Time Domain

In time domain, multiple time windows distributed at equal intervals are designed, and the integrated average value of the received signal in each time window is calculated to obtain the final response curve. The induced voltage after windowing procedure in time domain is obtained:

$$V_w\left(t\right)={\displaystyle {\int }_{-\infty }^{\infty }{V}_s\left(t\right)w_n\left(t-t_n\right)dt}$$

(7)

$V_s\left(t\right)$ denotes the induced voltage signal after stacking, $V_w\left(t\right)$ denotes the voltage signal after windowing, $t_n$ is the center of the window function, $w_n\left(t\right)$ is the rectangular window function, written as:

$$w_n\left(t\right)=\left\{\begin{array}{cc}\frac{1}{T_n}\hfill & \left|t\right|\le \frac{T_n}{2}\hfill \\ 0\hfill & \left|t\right|>\frac{T_n}{2}\hfill \end{array}\right.$$

(8)

where $T_n$ is the length of the rectangular window. Analyzed in frequency domain, the Fourier function of the rectangular window is equivalent to a low-pass filter, so the windowing process can eliminate out-of-band noise. Figure 7 shows the comparison of received voltage signals before and after windowing, indicating that the late SNR of the response signal has been improved a lot through the time window processing.

3.4. The Inversion of EM Data

Similar to the inversion of magnetic field data, the secondary field generated by the target is approximated by a three-dimensional magnetic dipole for the long-distance detection system. The difference is, the magnetic moment of the dipole is a constant in the inversion of magnetic data, while the magnetic moment of the dipole in EM detection is a time-dependent parameter, because the primary field is time-varying.

According to the theory proposed by L. R. Pasion [23,34], the magnetic moment $m_{EM}\left(t\right)$ has the form:

$$m_{EM}\left(t\right)=\frac{2\pi }{{\mu }_0}{A}^T{F}^{EM}\left(t\right)A\cdot b_{EM}^P=h_{EM}\left(\phi ,\theta ,L_1(t),L_2(t),L_3(t)\right)$$

(9)

where A is the Euler rotation tensor, determined by the azimuth $\phi$ and dip $\theta$ that are the rotation angle and the inclination angle between the local coordinate system where the target is located and the global coordinate system. $b_{EM}^P$ is the primary field at the location of the target, $F^{EM}\left(t\right)$ is the electromagnetic polarization tensor, characterized by three diagonal elements $\left(L_1\left(t\right),L_2\left(t\right),L_3\left(t\right)\right)$, $L_i\left(t\right)$ (i = 1,2,3) denotes the intrinsic response of the 3D magnetic dipole along the i-th direction, and $L_i\left(t\right)$ is estimated as [34]:

$$\frac{\partial L_i}{\partial t}=k_i{t}^{-{\beta }_i}\exp \left(-t/{\gamma }_i\right)$$

(10)

where $k_i$ is the magnitude of the polarization, ${\beta }_i$ defines the power law decay, ${\gamma }_i$ is the fundamental time constant. Once the three parameters $k_i,{\beta }_i,{\gamma }_i$ are determined, the intrinsic response curve of the magnetic dipole versus time can be estimated.

Combined (9) and (10), the function $h_{EM}$() is a unified representation for Equations (9) and (10), so the magnetic moment of dipole is described as:

$$m_{EM}\left(t\right)=h_{EM}\left(\phi ,\theta ,k_1,{\beta }_1,{\gamma }_1,k_2,{\beta }_2,{\gamma }_2,k_3,{\beta }_3,{\gamma }_3\right)$$

(11)

Based on the hypothesis of magnetic dipole, the EM feature vector y contains 14 unknown parameters, including 3 position parameters related to the location of target, 2 angle parameters related to posture of target and 9 parameters related to shape characteristics of target:

$$y=\left(x_0,y_0,z_0,\phi ,\theta ,k_1,{\beta }_1,{\gamma }_1,k_2,{\beta }_2,{\gamma }_2,k_3,{\beta }_3,{\gamma }_3\right)$$

(12)

The objective function of the electromagnetic inversion is established to minimize the misfit between the received voltage signal $d_{obsi}$ and the predicted the dipole response $g_i\left(y\right)$. The function $g_i$() in (13) is a unified representation for calculating the magnetic field generated by a dipole in Equation (1), the objective function of EM inversion is derived as follows:

$$\arg \underset{y}{\min }{\displaystyle \sum _{i=1}^N{\left(g_i\left(y\right)-d_{obsi}\right)}^2}$$

(13)

Levenberg–Marquardt (LM) algorithm is adopted to invert the feature vector y, which is a trust region method for solving the nonlinear problem [35]. LM is a local optimization strategy with the advantage of fast convergence.

Figure 8 is the flow chart of the joint inversion. It is worth emphasizing that the inversion results of the magnetic field is set as the initial value of the EM inversion, providing a priori information about the target location. The joint inversion method belongs to the “feature level” sensor fusion [27], the inversion result of the joint magnetic field and electromagnetic data is a fused feature vector, containing a three-dimensional position and two attitude angles of the target, and the shape feature of the target can be inferred from the intrinsic response curve of the magnetic dipole, realizing the classification of the target.

4. Experiments and Results

We have conducted two field tests to evaluate the performance of the joint detection system, and performed the cooperative inversion of the magnetic and electromagnetic data. In the first experiment, the joint detection system worked in full coverage survey mode. Both UAVMAG and TDEM-Cart were employed to survey the whole test area, and the mapping result of magnetic field data and electromagnetic data were compared. In the second experiment, the joint detection system worked in the cued survey mode. The whole test area was surveyed by the UAVMAG firstly. Then, a series of small-area fine surveys in the vicinity of the targets were carried out by the TDEM-Cart system, the cooperative inversion has been applied for the location and recognition of target.

4.1. Full Coverage Survey

4.1.1. Experimental Setup

This experiment was conducted in a suburb, Hebei Province, China. Three metal targets were pre-buried into the test area, as shown in Figure 9a. Among them, 1# and 2# targets are made of iron, 3# target is an aluminum triangle target.

The joint detection system worked in the full coverage survey mode. Firstly, the magnetic data was collected by the UAVMAG. The flight profiles have been planned in advance, the profiles are along the north–south direction, the height above the ground is 2 m, the interval is approximately 0.5 m, the error is less than 0.2 m affected by the weather. There are a total of 40 lateral lines, covering the entire test area. Then, the EM data was collected by the TDEM-cart system, and 9 lines with spacing 1 m were planned to cover the entire test area. Figure 9b,c show the actual UAVMAG flight trajectory and TDEM-cart test trajectory, the turning tracks are cut, leaving only the parallel profiles.

The magnetic data and EM data have been collected by the UAVMAG and TDEM-Cart system. Two data sets are synchronized based on the spatial location information measured by RTK GPS, and the positioning error is at the centimeter level. A unified cartesian coordinate system is established, and the origin is the take-off point of the drone, the x-axis points to the east, the y-axis points to the north, the z-axis points to up.

4.1.2. Results

The raw magnetic data measured by the magnetometer 1 and magnetometer 2 are given in Figure 10a,b, the magnetic anomaly field is extracted from the total magnetic field by the signal processing methods in Section 3.1, and Figure 10c demonstrates the mapping result of the magnetic anomaly signal over the entire test area. Compared with the magnetic map of raw data, the SNR of the processed signal has been improved a lot. Two magnetic anomalies are clearly observed from Figure 10c, which are in line with the preset positions of the 1# and 2# targets, whereas 3# target made of aluminum cannot be detected from the mapping result of magnetic data. The experimental results are consistent with the theoretical analysis, UAVMAG can only detect targets with high permeability.

The electromagnetic signal received by the TDEM-cart system decays with time. In order to map the secondary field, we obtained the spatial distribution of the voltage signal at the receiver through stacking and windowing methods described in Section 3.3. Figure 11a,b are the map of the response signal measured by No.1 receiving coil at t = 0.67 ms and t = 4.5 ms. Both early and late response signals show that there are three obvious abnormal signals in the test area, indicating that three metal targets may be buried. Furthermore, the spatial location of the abnormal signal is in agreement with the position of pre-buried targets. It is worth noting that the 3# target can be detected by the TDEM-Cart system, avoiding missed detection.

Therefore, when the joint detection system works in the full coverage survey mode, the combination of magnetic and electromagnetic signal is able to improve the reliability of the detection results and avoid missed inspections.

4.2. Cued Survey

4.2.1. Experimental Setup

This experiment was also conducted in a suburb, Hebei Province, China. Three ferromagnetic targets were pre-buried into the test area, marked as the 4#, 5#, and 6# target as shown in Figure 12. The responses of 4# target and 5# target are similar to UXO, which are used to simulate the real UXO, and 6# target is used to simulate the clutter. 4# target was buried at an angle of 45 degrees in the north–south direction, 5# and 6# targets were buried vertically, the location of the three targets were measured by the RTK GPS in advance. The position, azimuth and dip, shape feature of three targets are given in

Table 1. The basic parameters of three targets.

	4# Target	5# Target	6# Target
Location (m)	(6.17,33.96,−0.8)	(6.31,26.61,−1.0)	(3.99,22.15,−0.1)
Attitude angle $\left(\phi ,\theta \right)$	$\left(0^{°},{45}^{°}\right)$	$\left(0^{°},0^{°}\right)$	$\left(0^{°},0^{°}\right)$
Shape feature	rod-like	rod-like	Irregualr-shape

.

The joint detection system worked in the cued survey mode, a unified coordinate system was established, realizing the synchronization of magnetic and electromagnetic data in space.

4.2.2. Results

Firstly, the raw magnetic data was collected by the UAVMAG, and signal processing methods in Section 3.1 have been executed to extract the magnetic anomaly field related to the target. Figure 13a displays the magnetic anomaly field over the entire test area. The magnetic anomaly field generated by an isolated target is bipolar, including a positive anomaly field and a negative anomaly field. Three targets in the test area are observed from the mapping result. The magnetic anomaly signal produced by a single target has been manually separated, as shown in Figure 13b,d. Then, the DE algorithm proposed in Section 3.2 is employed to invert the feature vector x, and R-squared (R²) is used to measure the convergence of the inversion model. The position and magnetic moment of the three targets and model fit are listed in

Table 2. The inversion results from magnetic data.

	4# Target	5# Target	6# Target
Location (m)	(5.91,34.07,−0.86)	(6.35,26.41,−1.07)	(3.85,21.81,−0.06)
Magnetic moment (Am2)	(−0.016,2.14,−2.39)	(−0.194,−2.21,2.54)	(−0.001,−0.4,−0.05)
Model fit	90.81%	94.79%	86.18%

.

Next, potential targets have been delineated in the test area according to feature vector x. TDEM-Cart system was used to carry out finer surveys over the local areas near the targets. Three maps of received voltage at t = 0.67 ms are depicted in Figure 14.

Finally, EM inversion is performed based on electromagnetic data in three local areas, the initial value of EM inversion is set according to the feature vector x, the feature vector y of three targets retrieved from EM data and model fit are given in

Table 3. The inversion results from EM data.

	4# Target	5# Target	6# Target
Location (m)	(6.11,34.00,−0.83)	(6.31,26.55,−1.03)	(3.89,22.08,−0.08)
Attitude angle $\left(\phi ,\theta \right)$	$\left(0^{°},{40.21}^{°}\right)$	$\left(0^{°},{9.86}^{°}\right)$	$\left({8.53}^{°},{12.88}^{°}\right)$
Model fit	91.79%	93.90%	96.82%

, and the magnetic dipole eigen responses of three targets are demonstrated in Figure 15.

The feature vector x inverted from the magnetic data only provides the position and magnetic moment of ferrous targets, the posture, shape cannot be estimated. Comparing the position parameters in Table 1 and Table 2, the positioning error of the individual inversion of magnetic data is evaluated. The minimum positioning error in X is 0.04 m, the maximum is 0.26 m; the minimum positioning error in Y is 0.11 m, the maximum is 0.34 m; the minimum positioning error in Z is 0.04 m, the maximum is 0.07 m. In addition, based on the empirical analysis, it can be considered that the 3# target with the smallest magnetic moment might be a smaller target, which is consistent with the actual situation. However, the size of 1# and 2# target cannot be judged from the inverted magnetic moment.

The feature vector y derived from combined magnetic and EM data provides more information, including the position of the target, the azimuth and dip, and the magnetic dipole eigen response. Comparing the position parameters in Table 1 and Table 3, the positioning error of joint inversion in X, Y, and Z are all less than 0.1 m, verifying that the positioning accuracy estimated by the joint inversion is higher than that of individual inversion technique.

Comparing the angle parameters in Table 1 and Table 3, the inverted the attitude angle of the 4# and 5# targets are more accurate, error is less than 10°, whereas the calculation error of the 6# target is slightly larger. The reason is related to the assumption of the dipole model in the inversion theory, it’s difficult to use merely azimuth and dip to characterize its posture information for the irregularly shaped target.

Based on the dipole eigen response curve in Figure 15, the intrinsic responses of 4# target and 5# target have similar characteristics that the dipole response along the long principal axis is much larger than the dipole response in the other two directions, and decays more slowly over time, the dipole intrinsic curves in the direction of the two short principal axis basically coincide, indicating that 4# and 5# target may be rod-like and axisymmetric, which are in agreement with the characteristics of UXO. The intrinsic responses of 6# target in the three directions also attenuate with time, but its attenuation curves and amplitude are completely different, which conforms to the characteristics of irregularly shaped clutter.

5. Discussion

This paper introduces the joint remote detection system, which combines the UAV magnetic survey (UAVMAG) and the vehicle-mounted time domain electromagnetic system (TDEM-Cart), and the cooperative processing flow is designed, fusing the feature vectors inverted from magnetic data and electromagnetic data. Compared with the single detection technology [4,5,6,7,8,9], existing airborne and ground-based joint detection modes [24,25,26], the combination method proposed in this paper has significant advantages in terms of positioning accuracy and detection efficiency.

Compared with the individual magnetic detection, the joint detection system is capable of detecting both ferromagnetic and non-ferromagnetic metal targets, overcoming the shortcomings that only ferrous targets can be detected by the magnetic survey. In the first experiment, the joint detection system worked in full-coverage survey mode, and the experimental results confirm that the joint detection system is able to detect the target made of aluminum.

Compared with the single electromagnetic method, the joint detection system can improve the detection efficiency. When the joint detection system works in the cued survey mode, the rapid and complete scanning detection is firstly conducted by the UAVMAG, and potential abnormal regions are delineated according to the result of magnetic survey. Then, local fine surveys in the vicinity of the targets are performed by TDEM-Cart system, avoiding covering the entire test area and greatly shortening the working time.

For the cooperative interpretation, the inversion of magnetic data is based on the global optimization algorithm DE, the inversion of EM data fuses the feature vector inverted from magnetic data and adopts the local optimization algorithm LM. The proposed joint inversion process can not only improve the positioning result of the target, but also avoid falling into the local optimal solution. In the second experiment, the maximum positioning error of the magnetic field inversion reaches 34 cm, and the target positioning deviation after the joint inversion is less than 10 cm. In addition, the posture and shape characteristics of the target are inferred based on the two attitude angles and dipole intrinsic response over time, which are basically in line with the characteristics of the actual buried targets. The experimental results indicate that the accuracy of cooperative interpretation is higher, which is of great significance to the classification and recognition of targets. Nonetheless, it should be noted that the proposed joint processing algorithm is only suitable for the inversion of an isolated target.

Compared with the joint detection system in [23] with the combination of handheld magnetic survey and ground-based TDEM, the combined detection method proposed in this paper has higher detection efficiency than the ground magnetic survey, which is more suitable for large-scale survey and rapid monitoring. In [18], the UAV-mounted TDEM system has been developed, and the field experiments were carried out to detect UXO. The combination of the UAV-borne magnetic detection system and the UAV-borne electromagnetic system is bound to have better performance, which is the focus of future research.

6. Conclusions

In this paper, we concentrated on the remote detection of near-surface targets and interested in estimating the position, shape, and attitude information of the target. The main contents could be concluded that:

(1)

The joint detection system has been developed, including the drone-based magnetic detection system UAVMAG and the vehicle-mounted time domain electromagnetic detection system TDEM-Cart. The magnetic field data and electromagnetic data are synchronized through the spatial position measured by RTK GPS.

(2)

The processing flow consists of four parts, magnetic field signal processing, the inversion of magnetic data, EM data processing, and the inversion of EM data. The electromagnetic data inversion in the fourth step depends on the feature vector obtained from the inversion results in the second step, which belongs to the “feature level” sensor fusion.

(3)

In the first field test, the joint detection system worked in the full-coverage survey mode, both ferromagnetic and non-ferromagnetic metal targets can be detected.

(4)

In the second field test, the joint detection system worked in the cued survey mode, the position, shape, and attitude information of the target are derived based on the cooperative processing methods. The positioning accuracy is less than 10 cm and shows satisfactory performance in the recognition of target.