Continuous broadband lightning VHF mapping array using MUSIC algorithm

Atmospheric Research 231 (2020) 104647

Contents lists available at ScienceDirect

Atmospheric Research journal homepage: www.elsevier.com/locate/atmosres

Continuous broadband lightning VHF mapping array using MUSIC algorithm ⁎

T

⁎

Tao Wanga,b,1, Lihua Shia,1, , Shi Qiua,1, , Zheng Suna, Yantao Duana a b

National Key Laboratory on Electromagnetic Environmental Effects and Electro-Optical Engineering, Army Engineering University, Nanjing, China National Institute of Defense Technology Innovation, Academy of Military Sciences PLA China, Beijing, China

A R T I C LE I N FO

A B S T R A C T

Keywords: VHF MUSIC algorithm Time reversal Lightning localization Multi-source localization

In order to improve the mapping quality of lightning discharge processes, a mapping approach by applying multiple signal classification (MUSIC) algorithm to direction-of-arrival (DOA) estimation of lightning very high frequency (VHF) radiation sources is proposed (MUSIC-VHF method) and validated. Based on the principles of the eigen decomposition of the spatial covariance matrix, broadband MUSIC algorithm operating in lightning VHF ranges, on arbitrary planar array, is derived. To verify this method, numerical simulations were conducted for testing its performance in localization accuracy, spatial resolution and sidelobe suppression. Field experiments were designed to validate the proposed method, with an unmanned aerial vehicle (UAV) and a portable radiation source. The proposed method was further applied to map two classical triggered flashes. Temporal and spatial development of the discharge channels from the two flashes could be well depicted, including the upward positive leaders (UPLs), K-events, dart leaders and VHF radiations in subsequent return stroke (SRS). Additionally, comparisons among the proposed method and EMTR-VHF method, interferometer (INTF) technology are conducted. Comparison results show that the MUSIC-VHF method performs better in spatial resolution, sidelobe narrowband interference suppression. Moreover, the phenomena of multi-source radiations occurring simultaneously within a single observation window are reported, and their localization results are presented.

1. Introduction Lightning breakdown process emits electromagnetic radiation over a broad and continuous spectrum of radio frequencies, covering from ULF/VLF to VHF/UHF and higher. VHF radiation (30-300 MHz) is produced during smaller-scale breakdowns and has been employed to investigate lightning discharges in detail. In the past decades, time-ofarrival (TOA) techniques (Rison et al., 1999; Krehbiel et al., 2008; Thomas et al., 2004; Edens et al., 2012) and interferometer (Shao et al., 1995; Zhang et al., 2008; Mardiana and Kawasaki, 2000b; Yoshida et al., 2010; Yoshida et al., 2012) techniques have been widely used to map lightning radiations in two or three spatial dimensions, which helps to reveal the physical mechanism and gain insight into the characteristic parameters of lightning discharges. The basic principle for locating lightning radiation sources is to estimate the phase or timedifference of arrivals between two antennas. The first broadband interferometer was proposed by Shao et al. (1996), and it has been developed and improved quickly afterward

(Mardiana and Kawasaki, 2000a; Morimoto et al., 2005; Qiu et al., 2009a). Signals from pairs of antennas are coherently correlated to determine the phase or time difference of arrival of the two signals, and the angle cosines of the source direction are derived. Measurements along two orthogonal baselines are usually used to determine azimuth and elevation angles of the VHF radiations. Therefore, a triangle array consisting of three elements is generally employed to form an interferometer array in most applications. Since phase ambiguity problem affects the reliability of location, the phase difference technique has gradually been replaced by the time difference estimation technique which can avoid this problem. The broadband interferometry cross-correlation method is employed by Qiu et al. (2009b) to map the development of a triggered lightning process. Cao et al. (2012) and Sun et al. (2013, 2014, 2016) extended the application of the broadband cross-correlation method and proposed a general correlation time delay estimation algorithm based on direct correlation method and wavelet transform. Using a generalized crosscorrelation technique, Stock, 2014and Stock and Krehbiel, 2014

⁎

Corresponding authors. E-mail addresses: [email protected] (L. Shi), [email protected] (S. Qiu). 1 This work was supported in part by the National Science Foundation of China under Grant 41675002, Grant 41905002 and Grant 51977219. https://doi.org/10.1016/j.atmosres.2019.104647 Received 14 December 2018; Received in revised form 22 July 2019; Accepted 9 August 2019 Available online 13 August 2019 0169-8095/ © 2019 Elsevier B.V. All rights reserved.

Atmospheric Research 231 (2020) 104647

T. Wang, et al.

system. The position of each element is known as Am = (xm, ym), m = 1, 2, 3,∙∙∙, M. Note that, in an azimuth/elevation system, ϖ will consist of θ and φ which are the elevation and azimuth, respectively. Thus, the multiple signal classification approach begins with the following model for characterizing the received M frequency domain vector X(f) at each frequency point as

designed a continuous broadband digital interferometer. It is proved that continuous recording is more effective than sequentially triggered system in retrieving weak radiations. The continuous broadband digital interferometer technique has also been employed in (Rison et al., 2016; Zhang et al., 2017; Stock et al., 2017). As described above, traditional interferometer method for lightning location uses a three-element array and obtains the azimuth and elevation angles from two baselines. If one of them was noise-polluted the localization cannot be performed. Considering this problem, Stock et al. (2014) used a multiple baseline system to detect the VHF sources. In his research, the multiple antennas array was divided into a series of baselines composed of two-element array, and the direction of arrival (DOA) estimated from each baseline was averaged to derive a more accurate value. Recently, an electromagnetic time reversal (EMTR) technique, without estimating the phase or time-difference, begins to be applied to lightning localization. Mora et al. (2012) and Lugrin et al. (2014) first employed EMTR technique to locate the lightning return striking point. Wang et al. (2017, 2018) applied EMTR technique to lightning VHF localization and proved its performance in multi-source localization through field experiments. In (Wang et al., 2018), EMTR technique was further applied to multiple VHF sensors and continuous recording system. The performance of multi-antennae array versus the three-element array was also proved. EMTR technique uses the superposition of back-propagated waveform to locate the source. The advantage of this method is that it hardly needs any pre-processing of the detected waveforms. From the viewpoint of array signal processing, DOA estimation still has some other choices. Multiple Signal Classification (MUSIC) algorithm is one of the super-resolution spatial spectrum estimation algorithms (Schmidt, 1986). Audone and Margri (2001) applied MUSIC algorithm to finitedistance sources without omnidirectional radiation patterns to identify their positions and radiation characteristics. Chen and Lin (2011) used MUSIC algorithm to automatically compensate for the phase center shift during calibrating the free-space antenna factors of broadband electromagnetic compatibility antennas. While this method has been successfully used in many detection and localization applications (Ascione et al., 2013; Dell'Aversano et al., 2017), the use of it in lightning location has not been reported before. The small number of traditional lightning VHF locating antennas limits the use of advanced signal processing methods. With the use of multiple antennas, the resolution of DOA finding and the noise-suppression ability can be further improved. MUSIC algorithm is a special direction of arrival (DOA) estimation method based on the eigen decomposition of the sample covariance matrix. Utilizing the orthogonality between the signal subspace and the noise subspace, the eigenvalue of covariance matrix can be used to separate the signal from noise and a special spatial spectrum of array signals could be constructed. Due to its high resolution and strong sidelobe suppression ability, MUSIC algorithm might be a prospective tool for lightning localization based on multiple antennas array. In this paper, we first apply the MUSIC algorithm to lightning VHF localization. With the continuous broadband lightning VHF mapping system, MUSIC-VHF method is proposed to map the whole lightning discharging processes. The basic principles of MUSIC-VHF are described at the beginning of the paper. Its performances are analyzed, then compared with the EMTR-VHF method and interferometer (INTF) technology. Furthermore, with the proposed method, the phenomena of multi-source radiations during a triggered lightning are first reported in the end.

⎡ S1 (fk ) ⎤ ⎡ N1 (fk ) ⎤ ⎡ X1 (fk ) ⎤ ⎢ S2 (f ) ⎥ ⎢ N2 (f ) ⎥ ⎢ X2 (fk ) ⎥ k k ⎥+⎢ ⎥ = [ a (fk , ϖ1) a (fk , ϖ2) … a (fk , ϖQ )] ⎢ ⎢ ⎥ ⎢ ⋮ ⎥ ⎢ ⋮ ⎥ ⎢ ⋮ ⎥ ⎢ SQ (fk ) ⎥ ⎢ NM (fk ) ⎥ ⎢ XM (fk ) ⎥ ⎦ ⎣ ⎣ ⎦ ⎣ ⎦ (1a) or in matrix form

( (

) ⎤⎥ ⎥ ) ⎥⎥

⎡ exp −j2πf x1 cos θq cos φq + y1 cos θq sin φq k c ⎢ ⎢ x2 cos θq cos φq + y2 cos θq sin φq ⎢ exp −j2πfk c a (fk , ϖq) = ⎢ ⎢ ⋮ ⎢ x cos θq cos φq + yM cos θq sin φq ⎢ exp −j2πfk M c ⎣

(

)

⎥ ⎥ ⎥ ⎦

X (fk ) = A (fk , ϖ ) S (fk ) + N (fk )

(1b) (1c)

where fk represents the kth sampled frequency point, k = 1, 2, 3,···, K; the incident signals are represented in amplitude and phase by the complex quantities S(fk) = [S1(fk), ···, SQ(fk)]T; X(fk) = [X1(fk),···, XM(fk)]T represents the received M complex vector at each frequency point; “T” denotes transpose; a(fk, ϖq) depends on the mth element's position relative to the origin of the coordinate system, and its response to a signal (at fk) incident from the direction of the qth signal. The qth column of A is a “model” vector (Schmidt, 1986) a(fk, ϖq) responding to the qth signal. Knowing the model vector a(ϖq) is tantamount to knowing ϖq. The noise, whether “sensed” along with the signals or generated internal to the instrumentation, appears as the complex vector N(fk) = [N1(fk),···, NM(fk)]T. Note that the problem of solving for the directions of arrival of multiple incident wavefronts is actually to find the right a(f, ϖq) corresponding to the qth incident source within the range space of A. (In practice, a(f, ϖi) is obtained by meshing the φθ–plane into W directions, i = 1, 2, 3,···, W.) 2.1. The R matrix The solving of A starts from the M × M covariance matrix (R) of the measured data

R (fk ) = X (fk ) X H (fk )

(2)

where “H” denotes complex conjugate transpose. According to the principle of classical MUSIC algorithm (Schmidt, 1986), R can be partitioned into two orthogonal subspaces, namely the signal subspace and the noise subspace, under the basic assumption that the incident signals and the noise are uncorrelated. Then, the eigenvalue decomposition of covariance matrix can be carried out as

R (fk ) = U (fk ) ΣU H (fk ) = US (fk ) ΣS USH (fk ) + UN (fk ) ΣN UNH (fk )

(3)

where U is eigenvector matrix; Σ is a diagonal matrix composed of eigenvalues; US and UN are signal eigenvector matrix and noise eigenvector matrix, respectively; ΣS is a diagonal matrix consisting of large eigenvalues; ΣN is a diagonal matrix consisting of small eigenvalues (toward zero).

2. Broadband MUSIC algorithm 2.2. The spatial angular spectrum The lightning VHF signals received at the M elements array are linear combinations of noises with the Q incident waveforms from directions ϖq (q = 1, 2, 3,∙∙∙, Q) relative to the origin of the coordinate

According to the orthogonality between noise subspace and signal subspace, the noise eigenvectors are orthogonal to the incident signal 2

Atmospheric Research 231 (2020) 104647

T. Wang, et al.

mode vectors a(ϖq), as Eq. (4a) shows. Then, the estimation of a(ϖq) is equivalent to find the maximum of the spatial angular spectrum P(fk, ϖi) shown in (4b) by searching all the possible mode vectors a(ϖi).

aH (ϖq ) UN = 0

P (fk , ϖi ) =

Fig. 1(b) shows the frequency spectrum of the two signals. In this figure, the narrow band FM interference around 100 MHz is intensive. In this study, we directly remove the FM frequency components during postprocessing. Finally, the frequency ranges of 25 MHz~90 MHz and 110 MHz~150 MHz are selected for the localization algorithm. Fig. 1(c) shows the angular spectrum P(ϖi) distribution on the φθplane which is meshed with grid size of 0.1° × 0.1°. Through searching the location of the P(ϖi) peak in the plane, the direction to the impulsive radiation event will be obtained. In order to ensure both computational efficiency and angular resolution, two-step spatial searching method in (Wang et al., 2018) is also used in this paper, that is, the angular resolution can reach 0.01° in practical application. The P(ϖi) image in Fig. 1(c) also shows that the MUSIC-VHF method can suppress sidelobes effectively and form a sharp beam which indicates that the proposed method has high robustness.

(4a)

1 aH (fk , ϖi ) UNH (fk ) UN (fk ) a (fk , ϖi )

(4b)

For all the frequency points of the broadband signal, the spatial angular spectrum at each frequency point is properly combined to obtain the cost function (Wax et al., 1984).

P (ϖi ) =

1 K

K

∑ k=1

1 aH (fk , ϖi ) UNH (fk ) UN (fk ) a (fk , ϖi )

(5)

Then the beam pointing to the DOA is selected by searching the peak position of the following function on the φθ–plane.

ϖS = arg{max[P (ϖi )]}

2.4. Noise reduction

(6) Continuous recording and a heavy overlap used in the process of implementation can produce numerous solutions which are either contaminated by or entirely produced by noise. In order to filter those noise solutions, two metrics developed in (Wang et al., 2018), coherence ratio (CR) and energy ratio (ER), are also employed in this paper. Details of these two metrics can be found in (Wang et al., 2018). CR characterizes the consistent continuity of the object's spatial movement over a short time horizon. For an actual lightning VHF radiation, it appears to slowly drift over a small region during a short duration in the sky, resulting in large CR. However, for a noise solution, they tend to shift large distances across the sky randomly, resulting in small CR. ER, which is related to signal amplitude, characterizes the concentration of the P(ϖi) image. Solutions with lower ER are likely to be eliminated, for it suggests that the signal is seriously contaminated by or entirely noise.

where ϖs represents the DOA of the source. In summary, the steps of applying MUSIC to estimate the direction of lightning VHF radiations are: Step0: collect data, implement DFT, form R; Step1: calculate eigenstructure of R using (3); Step2: evaluate P(ϖi) versus ϖi using (5); Step3: pick peaks of P(ϖi) using (6). According to the above derivations, using aH(ϖq)UN = 0, the spatial spectrum constructed by MUSIC algorithm is similar to some kind of inversion spectrum. That is, in the source direction, P(ϖi)∝(1/0) → ∞, which suggests that MUSIC-VHF can produce sharp pencil beam shapes, resulting in strong sidelobe suppression and high spatial resolution. 2.3. Implementation The configuration of the continuous broadband lightning VHF mapping system used in this paper is the same with that described in (Wang et al., 2018). The lightning VHF radiations were recorded with 12-bit resolution at a sampling rate of 500 MHz continuously for 500 ms with 250 ms pretrigger time by an eight-channel LeCroy digital storage oscilloscope (DSO). The recorded signals were postprocessed utilizing MUISC-VHF method and applied to successively overlapping data segments. For the triggered lightning recorded in this paper, each segment has 256 samples (512 ns), and the successive window is displaced 64 samples (128 ns) relative to the preceding one. For a single window including an impulsive radiation event (also used in (Wang et al., 2018)) of the triggered lightning, the process of determining P(ϖi) distribution on the φθ–plane is depicted in Fig. 1. Fig. 1(a) shows two channel waveforms within a 256-sample window, with the vertical scale representing normalized amplitude. It is apparent that channel 6 (red) is advanced relative to channel 2 (blue).

3. Simulations Numerical simulations were conducted to investigate the performances of the MUSIC-VHF method in locating lightning VHF radiations, including localization accuracy under different background noise, spatial resolution, and sidelobe suppression. Comparisons with EMTRVHF method in (Wang et al., 2018), and INTF technology in (Stock, 2014; Zhang et al., 2017), are also included. Seven elements in the antenna array are configured in an “L” shape which is the same as the one employed in (Wang et al., 2018), as Fig. 7 shows. The waveform model is set to be a Gaussian pulse modulated oscillation of which the function is shown in (7).

s (t ) = A0 sin(2πf0 t ) e

‐4π

2

( t‐ττ ) 1 2

(7)

Fig. 1. A temporally isolated radiation event of the triggered lightning of section V, illustrating the processing technique used to determine P(ϖi) distributing on the φθ-plane. (a) Normalized waveforms from the two channels for a 256 (512 ns) sample window; (b) spectrum of the two signals; (c) distribution of normalized P(ϖi) for this isolated radiation event by the MUSIC-VHF method. 3

Atmospheric Research 231 (2020) 104647

T. Wang, et al.

Fig. 2. Robustness simulation results: (a) RMSE versus SNR; (b) ER versus SNR for MUSIC-VHF method.

EMTR-VHF when SNR ≥ -4 dB.

where f0 is the center frequency; τ1 and τ2 are the parameters that regulates the bandwidth of the signal; A0 is the amplitude of the signal. During the simulation, we set f0 = 100 MHz, τ1 = τ2 = 30 ns, and the frequency range is 0–200 MHz.

3.3. Narrow band interference suppression In a real scenario, FM interference (88 MHz~108 MHz) is mixed into the received signals. In order to quantify the effect of FM interference on the localization methods, simulations are performed as the ambient background noise is FM interference. And a parameter, named signal-interference ratio (SIR), is defined as the ratio of peak power of the interference and signal to rank the level of the interference, as shown in (9).

3.1. Accuracy analysis Assuming the noise is Gaussian noise with zero mean and variance σ . A far field broadband signal impinges on the array from ϖs = (120°, 45°). And white Gaussian noise with different SNR levels are added independently to each channel. After 10,000 independent simulations, the RMSE for MUSIC-VHF method was calculated from (8) and shown in Fig. 2(a). Meanwhile, the statistical mean result of ER versus SNR is shown in Fig. 2(b). 2

RMSE =

1 10000

SIR =

n=1

(9)

where PS represents the peak power of signal; PI represents the peak power of FM interference. The FM interference data used in the simulations was acquired from the observation site. The influence of FM interference on the localization results can be investigated by adjusting the SIR value. During the simulations, SIR increases from 0 to 1 at an interval of 0.01. The RMSE versus SIR, for MUSIC-VHF, EMTR-VHF and INTF methods are shown in Fig. 5. As Fig. 5 shows, the FM interference suppression of MUSIC-VHF method is superior to that of EMTR-VHF and INTF methods. If RMSE≤1° is required, the minimums of SIR for three methods should be 0.08, 0.24 and 0.31 respectively. The presence of narrow band interference causes the sidelobe level elevating in the distribution of P(ϖi), resulting in the decrease of localization accuracy. The advantage of the MUSIC-VHF method in antiinterference is attributed to its stronger sidelobe suppression (as mentioned earlier). For instance, when SIR = 0.3, the P(ϖi) imaging results of EMTR-VHF method and the proposed method are shown in Fig. 6. As Fig. 6 shows, MUSIC-VHF method performs much better in sidelobe suppression. Because MUSIC-VHF uses the zero point of the denominator in (5) to identify the DOA, its peak in the spatial angular spectrum is much sharper and higher than the peak obtained in EMTRVHF by direct accumulation of arrived signals.

10000

∑

PS PI

‖ϖn − ϖs ‖2 (8)

where ϖn is the nth estimation results at each SNR level. As shown in Fig. 2(a), for the white Gaussian noise case, the localization accuracy of MUSIC-VHF and EMTR-VHF methods is basically the same. For RMSE≤1°, the SNR of the two methods can be as low as -12 dB, which indicates high robustness of the algorithms. In contrast, the corresponding noise level for INTF is -4 dB. The INTF technology is implemented with antennas V1, V4 and V7 being used. Fig. 2(b) shows that ER increases with the increasing of SNR. An appropriate threshold for ER can be set to distinguish between actual lightning sources and noise events. For example, if RMSE≤1° is required (corresponding to SNR = -12 dB), the threshold for ER should be set at least to 0.18. 3.2. Spatial resolution Comparisons are also made to evaluate the spatial resolution of MUSIC-VHF and EMTR-VHF. Fig. 3 (a-b) show the localization results for source ϖs = (120°, 45°) by the two methods at SNR = 10 dB, respectively. In Fig. 3 (c), half peak width of main lobe is defined as the spatial resolution. For instance, the spatial resolutions of azimuth (SRA) for the two methods are 3.34° and 6.16°, respectively. As shown in Fig. 3, directive property of spatial spectrum resulted from MUSIC-VHF method is better than that by EMTR-VHF method. It suggests that MUSIC-VHF method can map the spatial structures of the lightning discharge channel in a finer way. In Section V, observation data is employed to evaluate this performance of the proposed method. Fig. 4 shows comparisons for spatial resolution of MUSIC-VHF and EMTR-VHF at different SNR levels. The spatial resolution curves for MUSIC-VHF and EMTR-VHF are very close to each other when SNR is from -12 dB to -6 dB. The SRA for MUSIC-VHF becomes better than

4. Experimental validation Employing an UAV and a portable radiation source, field experiments were designed and carried out to test the performance of MUSICVHF method with the “L” array mentioned earlier. 4.1. Instruments The configuration of the VHF array is shown in Fig. 7 (a). The incident signal is picked by seven VHF broadband omnidirectional 4

Atmospheric Research 231 (2020) 104647

T. Wang, et al.

Fig. 3. Spatial resolution comparison results at SNR = 10 dB: (a) by EMTR-VHF method; (b) by MUSIC-VHF method; (c) SRA comparison.

marked by red circles in Fig. 7 (b). At each point, a low elevation angle and a high elevation angle were tested. 4.2. Results and analysis At each test point, experiments were repeated three times independently. Taking one from each group of the tests, the received signals and their normalized P(ϖi) obtained by MUSIC-VHF method and EMTR-VHF method are shown in Fig. 8. Localization results by the two methods for the whole test data are shown in Table. 1. As shown in Table 1, the localization accuracy obtained by these two DOA estimation methods are generally accordant with each other. However, if we have a look of the distribution of P(ϖi), the advantage of MUSIC-VHF method becomes very clear. As Fig. 8 shows, MUISC-VHF method produces sharp main lobes with lower sidelobes and better spatial resolution. Note that the results in Fig. 8 and Table 1 are obtained after eliminating FM interference, as mentioned in section II(D). Taking TP-2 (78.4°, 12.3°) for an example, without removing the narrowband interference, the results by the proposed method and EMTR-VHF method are shown in Fig. 9. As Fig. 9 shows, when the FM interference is retained, the mapping result by EMTR-VHF method presents high sidelobe level which is about 65% of the main lobe, however, MUSIC-VHF still shows better sidelobe suppression ability. The conclusions derived from the field experiments and the simulations are in good agreement, which suggests that the proposed method, with stronger sidelobe suppression, higher spatial resolution and better anti-interference ability, possesses a promising prospect in accurate localization of lightning discharge processes.

Fig. 4. Spatial resolution comparison results at different SNR level.

Fig. 5. FM interference suppression results for MUSIC-VHF, EMTR-VHF and INTF methods.

5. Observations of triggered lightning and results

antennas with baseline of 9 m, and then recorded at a rate of 1.25GHz with 12-bit resolution by an eight channel LeCroy DSO. In Fig. 7 (b), a schematic diagram shows the position of the VHF array (black triangle) and the triggered lightning launch site (black pentagon) 1.1 km from the array. During the experiments, the UAV was used to carry a portable pulse source to a predetermined altitude and hovered in the air, meanwhile its azimuth and elevation were measured by a theodolite. The pulse source is connected to a radiating antenna in the form of coils, as shown in Fig. 7 (b). Two test points (TPs) in different orientations of the array were selected to in the test, TP-1 in the north and ~210 m away from the array, and TP-2 in the south and ~400 m from the array, which are

In this section, observation data for two triggered flashes are applied to validate the performance of MUSIC-VHF in practical applications. The discharges were triggered by classical artificial trigger technique. The observation system was installed at a site 1.1 km away from the rocket launcher, with 188.35° azimuth, as shown in Fig. 7, and its set was the same as (Wang et al., 2018). With direction finding methods, this continuous broadband lightning VHF mapping system can map lightning in two spatial dimensions with a time resolution of 0.128 μs. And the two metrics, CR and ER, mentioned in section II(D), are utilized to filter the noise solutions. The first flash case employed in this paper is the same as that in (Wang et al., 2018). The discharge was triggered without subsequent return stroke, called as slow discharges (Hubert and Laroche, 1984). After the flash being initiated, UPL developed upward from the rocket to propagate into the cloud, then an initial continuous current (ICC) 5

Atmospheric Research 231 (2020) 104647

T. Wang, et al.

Fig. 6. P(ϖi) imaging results by MUSIC-VHF and EMTR-VHF methods at SIR = 0.3: (a) by MUSIC-VHF method; (b) by EMTR-VHF method.

followed and lasted to the end. The second flash case is a classical negative triggered flash with obvious subsequent return strokes. Actually, the first case was triggered following the second case with an interval of five minutes. The amplitude statistics for the two triggered events are shown in Fig. 10. It can be found that the amount of large amplitude of the first case is smaller than that of the second case. The localization results presented for the two triggered events are also compared with the results by EMTR-VHF in (Wang et al., 2018), and INTF technology in (Stock, 2014; Zhang et al., 2017). Furthermore, multiple sources within the same time window are found in this triggered lightning. It is the first time that the phenomenon of multi-source radiation is discovered in the observation data.

28.78% of the number of initial solutions. However, according to (Wang et al., 2018), this percentage resulting from EMTR-VHF method is 26.99%. It suggests that MUSIC-VHF method produces more solutions passing CR threshold. The statistical results of ER of the proposed method are presented in Fig. 11 (c). Combining to the simulation results of ER in Fig. 2 (b), the threshold ER = 0.3 is selected to obtain a better mapping result. The final version for the results by MUSIC-VHF method is shown in Fig. 12 (a), in which the propagation paths become very clear after removing the noise points. All the data points are colored according to time, similarly hereinafter. For comparison purposes, the mapping result (Wang et al., 2018) by EMTR-VHF and INTF methods are shown in Fig. 12 (b) and (c). The rectangular box in each figure shows the number of remaining solutions. Similar to the above processing, the mapping results for the second lightning case by these three methods are shown in Fig. 12 (d), (e) and (f), respectively. Comparing the results in Fig. 12 (a), (b) and (c), it illustrates that: (i) the locatable VHF sources for the first lightning case by MUSIC-VHF method is the most among the three results; (ii) the details of the development of UPLs in this lightning are well depicted by the proposed method and EMTR-VHF method, but some of the discharge branches are not mapped clearly by INTF. Carefully comparing Fig. 12 (a) and (b), one can find that the MUSIC-VHF method can provide clearer imaging quality for some branch structures than EMTR-VHF method, for instance, branches L1 and L2 in Fig. 12. For the second flash case shown in Fig. 12 (d), (e) and (f), this triggered lightning can be contoured consistently by both the three methods. Nevertheless, the locatable VHF radiations by the proposed method are still the most among the mapping results.

5.1. Results and comparisons

y

29.5°

Taking the first lightning case for example, Fig. 11 (a) is the raw localization results prior to any metric filtering and shows the noisiness of the observations resulting from data recording being continuous. The statistical results of CR and ER are presented in Fig. 11 (b) and (c). Fig. 11 (a) presents all the solutions in the entire sky before any noise effects has been removed. Note that the mapping result for the first flash case is shown in the direction planar projection using cosα = cosθcosφ, and cosβ = cosθsinφ, similarly hereinafter. The progresses of reducing noise solutions are the same as the two phases depicted in (Wang et al., 2018). In phase 1, the metric CR is used to filter false solutions produced by noise. And the threshold of CR is also set as 0.3. In phase 2, the metric ER is employed subjectively to filter the remaining noise solutions. As Fig. 11 (b) shows, the number of solutions with CR ≥ 0.3 take up

z

8.35°

(a)

(b)

Fig. 7. Schematic diagram of experimental setup: (a) the antenna array; (b) the radiation source and points. 6

Atmospheric Research 231 (2020) 104647

T. Wang, et al.

Fig. 8. UAV experimental results by MUSIC-VHF and EMTR-VHF: (a), (d), (g), and (j) are the received signals from different radiation sources; (b), (e), (h), and (k) are the mapping results by the proposed method; (b), (e), (h), and (k) are the mapping results by EMTR-VHF method.

Actually, the above conclusions are consistent with the simulation results in Fig. 2 (a), Fig. 5, Fig. 6 and the experimental results in Fig. 9. It indicates that INTF technology is more susceptible to the environment noise, especially for the case of weak radiations. As for MSUICVHF and EMTR-VHF methods, it can be inferred that the advantage of MUISC-VHF in locating faint radiations benefits from its stronger

Taking further comparisons of performances of the three methods in locating faint sources, for each lightning case, Fig. 12 (g) and (h) show the statistical results of the received signal amplitude corresponding to each solution by MUSIC-VHF, EMTR-VHF and INTF methods, respectively. Both figures illustrate that the number of weak sources located by MSUIC-VHF method is the most among the methods. 7

Atmospheric Research 231 (2020) 104647

T. Wang, et al.

Table 1 Statistical results of UAV experiments by MUSIC-VHF and EMTR-VHF methods. Measurements (°)

TP1

TP2

(241.5, 14.7) (241.7, 14.5) (242.3, 15.3) (241.5, 30.5) (241.7, 31.1) (244.1, 28.5) (77.6, 8.0) (77.7, 8.1) (77.6, 7.9) (78.4, 12.3) (78.4, 12.4) (78.4, 12.4)

Localization results (°) MUSIC-VHF

Error

EMTR-VHF

(241.05, 14.2) (241.22, 14.33) (241.23, 15.7) (240.85, 30.24) (241.30, 30.57) (243.83, 28.83) (77.37, 8.10) (77.14, 7.46) (77.45, 7.40) (78.04, 11.67) (77.31, 11.74) (77.94, 12.02)

(0.45, 0.50) (0.48, 0.17) (1.07, 0.40) (0.65, 0.26) (0.4, 0.53) (0.27, 0.33) (0.23, 0.10) (0.56, 0.64) (0.15, 0.50) (0.36, 0.63) (1.09, 0.66) (0.46, 0.38)

(240.54, 14.10) (240.95, 13.69) (241.38, 16.1) (240.82, 30.24) (241.29, 30.41) (243.54, 29.11) (77.41, 7.12) (77.17, 7.01) (77.30, 7.10) (78.12, 11.62) (78.06, 11.88) (77.36, 12.09)

Error (0.96, (0.75, (0.92, (0.68, (0.41, (0.56, (0.19, (0.53, (0.30, (0.28, (0.34, (1.04,

0.60) 0.81) 0.80) 0.26) 0.69) 0.61) 0.88) 1.10) 0.80) 0.68) 0.52) 0.31)

Fig. 9. TP-2 (78.4°, 12.3°) localization results by MUSIC-VHF and EMTR-VHF methods without removing FM interference.

After about 6.9 ms later than the initiation, UPL developed away from the flash origin into multiple directions, the copious positive breakdown events were detected (Fig. 13 (i)). These faint radiations provide a good case to verify the performances of different methods in resisting external interference. As shown in Fig. 13 (i), in the duration of UPL developing upward, a succession of fast retrograde negative leaders (“K”-events, also referred to as “recoil” leaders) were also detected. The K-events developed back along inferred positive leader paths and were shown as vertical lines in the elevation versus time panel. One of the typical K-events will be analyzed in B and the results are shown in Fig. 15. As for the second lightning case, immediately following initiation, UPL developed upward. According to the mapping results, there are three return strokes events which were initiated by three dart leaders being contained in the observed data. Fig. 14 presents the three dart leaders' development versus time, using the mapping results by MSUICVHF method. Dart leaders propagate rapidly and produce copious VHF radiations, therefore their results reproduced by the three methods are similar with each other. However, when careful comparisons are taken, new information emerges after Fig. 14 is stretched. Some VHF radiations are detected in the subsequent return stroke (SRS) channel, that is considered to be in small amount and much weak. The localization results of VHF radiations in SRS by EMTR-VHF and INTF methods are also listed in Fig. 14. It illustrates that the result of EMTR-VHF is basically the same as that of the proposed method except that a few extra points at the upper right corner are captured by the latter. Whereas, VHF radiation sources in this breakdown discharge channel detected by INTF technology are far less and non-continuous.

Fig. 10. Amplitude statistics for the two triggered events.

narrowband interference suppression. For observation data, although the FM interference band with big amplitude is eliminated, other unknown or intermittent narrowband interference frequencies may still be hidden in the acquired array signal. Further, with respect to the first case, Fig. 13 shows the flash development in views of elevation/azimuth angle versus time. Different stages of this discharging processes are zoomed in Fig. 13 (iii)-(v) respectively. The initial 244 ms of activity depicts breakdown discharge events during the rocket ascending, and its rising angle is ~1.34° (1.36° for EMTR-VHF method) as shown in Fig. 13 (iii). The dashed red circles in Fig. 13 (i) and (ii) indicate where this triggered lightning was initiated. From Fig. 13 (ii) and (iii), it can be found that the flash was initiated 244 ms later at ϖ = (188.25°, 30.01°) ((188.2°, 30.11°) for EMTR-VHF method) and lasted to 430 ms. As Fig. 13 (i) and (iv) present, the processes of wire fusing are also revealed clearly. 8

Atmospheric Research 231 (2020) 104647

T. Wang, et al.

Fig. 11. Raw observation data and CR, ER results for the first flash case: (a) observations prior to any metric filtering; (b) CR results; (c) ER results.

Fig. 12. Mapping results by MUSIC-VHF, EMTR-VHF and INTF methods: (a), (b) and (c) for the first case; (d), (e) and (f) for the second case; (g) and (h) are statistical results of the receiving signal amplitude for each solution from the two cases respectively.

9

Atmospheric Research 231 (2020) 104647

T. Wang, et al.

Fig. 13. Triggered lightning development versus time: (i) elevation angle versus time; (ii) azimuth angle versus time; (iii)-(v) zooms of the stages of rocket rising and wire fusing.

Fig. 14. Dart leaders and VHF radiations in SRS mapping results by MUSIC-VHF, EMTR-VHF and INTF methods.

selected for comparison, for VHF radiations produced by negative breakdown during lightning discharge are generally stronger than that from positive leaders. The selected K-event is mapped at the same scale by the two methods in Fig. 15 (a) and (b), respectively, in which the white substrates are the contours of the triggered lightning. As the numbers in white boxes show, the number of radiations located by the two methods

5.2. Comparisons of spatial resolution According to the simulation results in Fig. 4 and field experimental results in Fig. 8, the spatial resolution of the proposed method is better than EMTR-VHF method in the case of high SNR. In order to further validate this conclusion, a section of observation data from the triggered lightning is selected for comparisons. One of the K-events is 10

Atmospheric Research 231 (2020) 104647

T. Wang, et al.

Fig. 15. Mapping results of the selected K-event: (a) by MUSIC-VHF method; (b) by EMTR-VHF method.

Fig. 16. Cases of multi-source radiations and their localization results: (a), (d) and (g) are data waveforms; (b), (e) and (h) are corresponding results by MUSIC-VHF method; (c), (f) and (i) are corresponding results by EMTR-VHF method. S1 and S2 represents two sources, respectively.

11

Atmospheric Research 231 (2020) 104647

T. Wang, et al.

Fig. 17. Mapping three cases of multi-source radiations in Fig. 16 in the UPL channels: (a)-(c) corresponding to (a), (d) and (g) in Fig. 16.

continuous recording mode with 500 MHz sampling rate. Noise filtering technology developed in (Wang et al., 2018), combining metrics CR and ER, is also applied to filter the false solutions. Comparisons among MUSIC-VHF method, the previous EMTR-VHF method and interferometer (INTF) technology are also carried out. Simulations and experimental researches were carried out to analyze the performances of the proposed method in localization accuracy, spatial resolution, sidelobe suppression and narrowband interference suppression. The statistical simulation results show that: (i) the localization accuracy of MUSIC-VHF and EMTR-VHF methods is basically consistent with each other in the case of white Gaussian noise, and both of them perform better than INTF; (ii) the spatial resolution of the proposed method is better than that of EMTR-VHF method in the case of high SNR; (iii) MUSIC-VHF method performs better than EMTR-VHF and INTF methods in suppressing sidelobe and resisting narrowband interference. Field experiments with an UAV and a portable radiation source were conducted. The robustness and validity of the proposed method are tested and verified in practical situations. Moreover, the proposed method is applied to observation data of a classical triggered lightning. Compared to EMTR-VHF and INTF methods, MUSIC-VHF method is able to reproduce the development of lightning at VHF in exceptional detail. Much more copious VHF radiations from positive leaders can be located. And some branch structures with relatively high amplitude are mapped with better fineness and convergence. This advantage is mainly due to the low side-lobe and sharp main lobe in DOA estimation (Fig. 3) by using MUSIC method. Furthermore, the phenomena of multi-source radiations occurring simultaneously within a single window are found in the observation data. With the proposed method, the two sources contained in the recording waveforms are recognized clearly. It indicates that the proposed method holds promise in being capable of observing the possible bidirectional leader development during lightning initiations. It needs to be noticed that, in selection of the localization method for lightning VHF radiation sources, one needs to make a compromise between the shortcomings and the advantages of the candidate methods. The INTF algorithm implementing with three antennas and the EMTR-VHF method are simpler and easier than the MUSIC-VHF method which requires eigendecomposition of covariance matrix. In addition, more antennas adopted in the array configuration will also increase the complexity of the system.

is basically equivalent. And the two imaging results are basically consistent in shape. It suggests that both the two methods work well in reproducing lightning discharge processes. However, after careful comparisons, the result by MUSIC-VHF method shows better fineness and convergence. 5.3. Multiple sources localization Lightning VHF multi-source localization has received attention and mentioned in some studies (Stock and Krehbiel, 2014; Stock et al., 2014; Wang et al., 2017; Wang et al., 2018). However, the results of lightning VHF multi-source localization are rarely reported, especially the multi-source radiation phenomenon in observation data. Additionally, multiple sources arriving at the observation system simultaneously, that are close to each other, cannot be distinguished by traditional direction-finding methods, including VHF interferometry and TOA techniques. As the lightning mapping result shown in Fig. 12 (b), the UPL in this flash developed in multiple branches. Some multiple simultaneous breakdown events occurred during the development of the triggered lightning. Fig. 16 (a), (d) and (g) present three raw waveforms from two simultaneous sources which arrive at the VHF array so closely that they could not even be identified artificially in a single window, especially for (a) and (d). The three events occurred at 252.316 ms, 258.628 ms and 316.618 ms, respectively. The best way to locate the two sources is to use multi-source localization methods. The second and third columns of Fig. 16 shows the dual-source localization results by MUSIC-VHF and EMTR-VHF methods, respectively. It illustrates that the two sources can be clearly recognized by the two methods, but the sidelobes in the results by EMTR-VHF method are stronger than that of MUSIC-VHF method. In Fig. 17, distributions of the three dual-source cases in the flash channels are presented. In this figure, the two sources (red dots) are connected by dash lines and the substrate is the contours of UPLs. To the best of our knowledge, it is the first time that the phenomena of multi-source radiations in observation data are reported and located. Although only a few isolated radiation events were found, it suggests that the proposed method hold promise in being capable of locating multiple branches occurring simultaneously. 6. Conclusions and summary

Acknowledgments

In this paper, a DOA estimation method by combining MUSIC algorithm with VHF observation is proposed, and the algorithm of MUSIC-VHF method is derived. Continuous broadband lightning VHF mapping system in (Wang et al., 2018) is still employed as an instrument. An antenna array consists of seven omnidirectional antennas with 9 m baseline in “L” shape. An eight-channel LeCroy DSO works in

This work is supported by the National Natural Science Foundation of China (Grant No.41675002; Grant No.41905002; Grant No.51977219). The authors would like to thank all the Jiangsu lightning observation team (SLOT) members taking part in the field 12

Atmospheric Research 231 (2020) 104647

T. Wang, et al.

lightning observation campaign.

10.1038/ncomms10721. Art. no. 10721. Schmidt, R.O., 1986. Multiple emitter location and signal parameter estimation. IEEE Trans. Antennas Propag. AP-34 (3), 276–280. Shao, X.M., Krehbiel, P.R., Thomas, J.R., Rison, W., 1995. Radio interferometric observations of cloud-to-ground lightning phenomena in Florida. J. Geophys. Res. 100 (D2), 2749–2783. Shao, X.M., Holden, D.N., Rhodes, C.T., 1996. Broad band radio interferometry for lightning observations. Geophys. Res. Lett. 23 (15), 1917–1920. Stock, M.G., 2014. Broadband Interferometry of Lightning. Ph.D. Dissertation. Sci. Phys., Langmuir Lab, New Mexico Inst. Mining Technol, Socorro, NM, USA. Stock, M.G., Krehbiel, P.R., 2014. Multiple baseline lightning interferometry–improving the detection of low amplitude VHF sources. In: Proc. Int. Conf. Lightning Protection (ICLP), Shanghai, China, pp. 293–300. Stock, M.G., Akita, M., Krehbiel, P.R., et al., 2014. Continuous broadband digital interferometry of lightning using a generalized cross-correlation algorithm. J. Geophys. Res. Atmos. 119 (6), 3134–3165. Stock, M.G., Krehbiel, P.R., Lapierre, J., Wu, T., Stanley, M.A., Edens, H.E., 2017. Fast positive breakdown in lightning. J. Geophys. Res. Atmos. 122 (15), 8135–8152. https://doi.org/10.1002/2016JD025909. Sun, Z., Qie, X.S., Liu, M.Y., et al., 2013. Lightning VHF radiation location system based on short-baseline TDOA technique—Validation in rocket-triggered lightning. Atmos. Res. 129–130, 58–66. Sun, Z., Qie, X.S., Jiang, R.B., et al., 2014. Characteristics of a rocket-triggered lightning flash with large stroke number and the associated leader propagation. J. Geophys. Res. 119 (23), 13388–13399. https://doi.org/10.1002/2014JD022100. Sun, Z., Qie, X.S., Liu, M.Y., et al., 2016. Characteristics of a negative lightning with multiple-ground terminations observed by a VHF lightning location system. J. Geophys. Res. Atmos. 121 (1), 413–426. Thomas, J.R., Krehbiel, P.R., Rison, W., et al., 2004. Accuracy of the lightning mapping array. J. Geophys. Res. 109 (D14), D14207. https://doi.org/10.1029/ 2004JD004549. Wang, T., Qiu, S., Shi, L.H., Li, Y., 2017. Broadband VHF localization of lightning radiation sources by EMTR. IEEE Trans. Electromagn. Compat. 59 (6), 1949–1957. Wang, T., Shi, L.H., Qiu, S., et al., 2018. Multiple-antennae observation and EMTR processing of lightning VHF radiations. IEEE Access 6, 26558–26566. Wax, M., Shan, T.-J., Kailath, T., 1984. Spatio-temporal spectral analysis by eigenstructure methods. IEEE Transactions on Acoustics ASSP-2 (4), 817–827. Yoshida, S., Biagi, C.J., Rakov, V.A., et al., 2010. Three-dimensional imaging of upward positive leaders in triggered lightning using VHF broadband digital interferometers. Geophys. Res. Lett. 37 (5), L05805. https://doi.org/10.1029/2009GL042065. Yoshida, S., Biagi, C.J., Rakov, V.A., et al., 2012. The initial stage processes of rocket-andwire triggered lightning as observed by VHF interferometry. J. Geophys. Res. Atmos. 117https://doi.org/10.1029/2012JD017657. no. D9, p. D09119. Zhang, G., Zhao, Y., Qie, X., Zhang, T., Wang, Y., Chen, C., 2008. Observation and study on the whole process of cloud-to-ground lightning using narrowband radio interferometer. Sci. China Ser. D Earth Sci. 55 (5), 694–708. Zhang, Y., Krehbiel, P.R., Zhang, Y.J., et al., 2017. Observations of the initial stage of a rocket-and wire-triggered lightning discharge. Geophys. Res. Lett. 44 (9), 4332–4340. https://doi.org/10.1002/2017GL072843.

References Ascione, M., Buonanno, A., D’Urso, M., et al., 2013. A new measurement method based on music algorithm for through-the-wall detection of life signs. IEEE Trans. Instrum. Meas. 62 (1), 13–26. Audone, B., Margri, M.B., 2001. The use of MUSIC algorithm to characterize emissive sources. IEEE Trans. Electromagn. Compat. 43 (4), 688–693. Cao, D., Qie, X.S., Duan, S., et al., 2012. Lightning discharge process based on shortbaseline lightning VHF radiation source locating system. Acta Phys. Sin. 61 (6), 069202. Chen, H.F., Lin, K.-H., 2011. An improved method for free-space antenna-factor measurement by using the MUSIC algorithm. IEEE Trans. Electromagn. Compat. 53 (2), 274–282. Dell’Aversano, A., Natale, A., Buonanno, A., et al., 2017. Through the wall breathing detection by means of a doppler radar and MUSIC algorithm. IEEE Sensors Letters 1 (3), 3500904. Edens, H.E., Eack, K.B., Eastvedt, E.M., et al., 2012. VHF lightning mapping observations of a triggered lightning flash. Geophys. Res. Lett. 39 (19), L19807. https://doi.org/ 10.1029/2012GL053666. Hubert, P., Laroche, P., 1984. Triggered lightning in New Mexico. J. Geophys. Res. 89 (D2), 2511–2521. Krehbiel, P.R., Riousset, J.A., Pasko, V.P., Thomas, R.J., 2008. Upward electrical discharges from thunderstorms. Nat. Geosci. 1 (4), 233–237. Lugrin, G., Parra, N.M., Rachidi, F., Rubinstein, M., Diendorfer, G., 2014. On the location of lightning discharges using time reversal of electromagnetic fields. IEEE Trans. Electromagn. Compat. 56 (1), 149–158. Mardiana, R., Kawasaki, Z., 2000a. Dependency of VHF broad band lightning source mapping on Fourier spectra. Geophys. Res. Lett. 27 (18), 2917–2920. Mardiana, R., Kawasaki, Z., 2000b. Broadband radio interferometer utilizing a sequential triggering technique for locating fast-moving electromagnetic sources emitted from lightning. IEEE Trans. Instrum. Meas. 49 (2), 376–381. Mora, N., Rachidi, F., Rubinstein, M., 2012. Application of the time reversal of electromagnetic fields to locate lightning discharges. Atmos. Res. 117, 78–85. Morimoto, T., Kawasaki, Z., Ushio, T., 2005. Lightning observations and consideration of positive charge distribution inside thunderclouds using VHF broadband digital interferometry. Atmos. Res. 76 (1), 445–454. Qiu, S., Zhou, B.H., Shi, L.H., Dong, W.S., Zhang, Y.J., Gao, T.C., 2009a. An improved method for broadband interferometric lightning location using wavelet transforms. J. Geophys. Res. 114 (D18), D18211. https://doi.org/10.1029/2008JD011655. Qiu, S., Yang, B., Dong, W.S., Gao, T.C., 2009b. Application of correlation time delay estimation in broadband interferometry for lightning detection (in Chinese). Scientia Meteorologica Sinica 29 (1), 92–96. Rison, W., Thomas, R.J., Krehbiel, P.R., Hamlin, T., Harlin, J., 1999. A GPS based threedimensional lightning mapping system initial observation in central New Mexico. Geophys. Res. Lett. 26 (23), 3573–3576. Rison, W., Krehbiel, P.R., Stock, M.G., et al., 2016. Observations of narrow bipolar events reveal how lightning is initiated in thunderstorms. Nat. Commun. 7https://doi.org/

13