Detection, Identification, and Direction of Arrival Estimation of Drone FHSS Signals with Uniform Linear Antenna Array

Safety, security, and privacy are three critical concerns affiliated with the use of drones in everyday life. Considering their ever-shrinking sizes and capabilities, being aware of drone activities in the vicinity becomes an important surveillance item. Therefore, keeping track of drones and preferably their controllers should be included into the already-existing security measures. In this study, a frequency hopping spread spectrum (FHSS) type drone controller signal detection and emitter direction finding framework is proposed to achieve aforementioned goals. Since drone communications signals generally coexist with other FHSS signals in 2.4 GHz industrial, scientific, and medical (ISM) band, first, a method based on cyclostationarity analysis is proposed to distinguish the drone radio controller signals from other signals utilizing 2.4 GHz ISM band. Then, a variant of short-term Fourier transform is introduced to extract the parameters of detected drone remote controller signals. The correct hopping signals are then aggregated based on the clustered parameters to obtain their combined baseband equivalent signal. Furthermore, the resampling process is applied to reduce the unrelated samples in the spectrum and represent the spectrum with the reconstructed signal, which has a much lower bandwidth than the spread bandwidth. Finally, two different multiple signal classification algorithms are utilized to estimate the direction of the drone controller relative to the receiving system. In order to validate the overall performance of the proposed method, the introduced framework is implemented on hardware platforms and tested under real-world conditions. A uniform linear antenna array is utilized to capture over-the-air signals in hilly terrain suburban environments by considering both line–of-sight and non–line–of-sight cases. Direction estimation performance is presented in a comparative manner and relevant discussions are provided.

other hand, these capabilities pose several concerns regarding safety, security, and privacy. Among these concerns, presence of drones is considered to be a threat, especially in the vicinity or surrounding airspace of critical zones [2]. In addition to being a direct physical threat to vehicles, infrastructures, and facilities, drones could be an auxiliary source of threat especially when they are employed to jam, interfere, or totally block communication links leading to malfunction in critical operations and/or services. It is evident that detection of any potential or actual intrusion by cooperative and non-cooperative drones and their pilots is required to protect critical zones, vehicles, operations, and infrastructures. Therefore, in this study, a method, which identifies the frequency-hopping spread spectrum (FHSS) type drone radio controller (RC) and estimates the angle of arrival (AoA) of the pilot of the drone, is presented. The proposed method is successfully implemented on hardware and tested with field measurements. Results and discussions are provided along with future directions.

A. RELATED WORK
Drone detection in a given area is an important task. One way of drone detection is to see if there are any drone controlling signals in a given environment. The first step to detect RC/pilot of a drone is to sense the presence of the control and/or communication signals exchanged between the drone itself and its RC. Majority of signals transmitted over the drone-RC link are of FHSS type since spread spectrum techniques provide resilience to interference, enhance security, and enable networking operations, as well as to provide covert transmission [3]. Hence, detection of a spread spectrum signal is one of the major challenges in locating the pilot of the drone. Even though drones exploit FHSS techniques, one should note that they operate in civilian airspaces by occupying industrial, scientific, and medical (ISM) bands shared by many other telecommunication links which employ FHSS as well. In this respect, there are various algorithms developed for signal detection over the years such as autocorrelation-based, wavelet-based, and timefrequency analysis based algorithms. Moreover, compressive sampling, which is recently developed with advanced digital signal processing techniques, is also adopted for FHSS signals detection [4,5]. One of the most well-known signal detection analysis is cyclostationarity which is used to detect the periodic signals [6]. In this regard, we also utilized the cyclostationarity analysis in order to distinguish signals. This poses a second challenge since the RC/pilot locating process should be able to discriminate and focus on drone-RC link signals rather than other non-threatening FHSS signals. Once a potential drone-RC/pilot link is detected, the next step is to determine the direction of the incoming signals. This immediately implies the employment of radio-direction finding (DF) techniques. Here, it is important to note that radio-DF has a rich history and a very well-established literature on wide variety of its application scenarios. Here, we will briefly visit some of the widely used approaches. DF techniques can be grouped under the following two categories: switched beam system (SBS) and adaptive array system (AAS) [7]. SBS method uses overlapping beams to scan the azimuth plane. The AoA is then determined by a search carried out across all over the candidate beam positions. Output is chosen to be the position yielding a maximum value for a defined cost function. It is important to state here that the SBS method benefits from the following two key points: firstly, it operates with a mechanically agile directional single antenna, and secondly, it does not require any baseband signal processing. Of course, these key points come at the expense of several shortcomings. For instance, the accuracy of SBS method relies heavily on signal-to-noise ratio (SNR). This points out that detection range affects the performance significantly. Furthermore, statistical nature of wireless propagation generally implies poor SNR regimes for links especially within regulated but unlicensed bands where strict power limitations are imposed on transmitters. Hence, optimal performance could only be achieved within the time slots when SNR is above a certain threshold. Therefore, the overall performance could be improved by prolonging the observation intervals, which is against the nature of the problem of interest. Such shortcomings necessitate more robust and responsive algorithms.
In contrast to SBS, AAS methods take the advantages of smart antennas to steer the main beam in any desired direction and establish a continuous tracking. Beam steering is achieved by combining the weights of each antenna array element so that the maximum power is transferred spatially for a desired location or direction. Equipped with beam steering/tracking capability, AAS could be used in DF with the following two groups of techniques: conventional and subspace-based. Conventional methods form a power spectrum in such a way that look-direction gain is forced to be at unity gain while those of all other directions are minimized [8]. This way, a peak search could be carried out across the spectrum for the candidate angles/directions for the signal source of interest. Here, it is noteworthy to indicate that conventional techniques could operate in low SNR regimes to some extent. However, their computational complexities are higher as compared to those of AAS. On the other hand, subspace-based techniques focus on separating the signal subspace orthogonal to the noise subspace. With the advances in digital signal processing techniques and systems, subspace-based methods could be implemented and deployed rapidly [9,10]. Among the various directionfinding algorithms, which have a different methodology, the correlative interferometer needs a calibration process to fill a phase offset table and requires an odd number of antennas. In addition to this, pseudo-Doppler has a resolution problem. Regarding expectation-maximization, while the result has a prominent accuracy but has difficult implementation problems, and also it requires prior information. Contrary, Multiple Signal Classification (MUSIC) algorithm gives results instantaneously, and also performs better under dynamically changing environments. Furthermore, MUSIC algorithm is robust under the multipath and mobility scenarios [11,12].
The majority of the DF algorithms present in the literature focus on a non-hopping carrier and operate on a narrowband spectrum [13,14]. When FHSS signals are considered, majority of the DF techniques could not be employed directly. The same reasoning holds for wideband signals as well. However, there are some efforts in the literature to lay down the implications of frequency changes in statistical characteristics of DF [15]. For instance, peculiar to drone RC localization, utilization of mechanically-agile directional antenna is proposed in [16]. However, [16] considers signals that are not of FHSS-type RC and thus the method's performance is not satisfactory for such signals. For the FHSStype drone RC, [17] proposed Anderson-Darling goodnessof-fit based passive signal detection and DF system. In [17], regardless of what type of a signal is present, it is not possible to distinguish the drone RC and hence the false alarm rate of the DF results can be high. There exists some non-linear methods in the literature as well. For instance in [18], a convolutional neural network (CNN) model is proposed to localize drone controllers where location information is obtained by a method that could be considered as radio frequency (RF) fingerprinting. Note that performance of this approach inherently relies on the data set to be used for training, which prevents it from being generalized. Considering the fact that RF environment is dynamic, the CNN is in need of retraining for every single change in any environment as well as for different type of drone controllers.

B. CONTRIBUTIONS
Discussions carried out up until this point reveal that a comprehensive approach, depending on both theoretical and real-world measurement data validation, is required for DF of RC signal sources since drone controllers manifest several interesting features simultaneously. Thus, the contributions of this study that exploit these features are: • To the best of authors' knowledge, this work provides one of the first FHSS type drone RC DF framework and measurement results for the real-world measurement campaign under line-of-sight (LOS) and non-line-ofsight (NLOS) scenarios.
• First, a method based on cyclostationarity analysis is applied to distinguish the drone RC signal from others. Furthermore, short time Fourier transform (STFT)based blind signal detection and resampling processes are employed to convert the captured wideband signal into a narrowband signal. Two of MUSIC variants are implemented at the final stage over the processed data to estimate the direction of arrival for drone RC.

•
The performance of the proposed method is empirically assessed and it is shown that the resultant framework can find the direction of arrival of the drone RC signal with an average of 1.39 degrees of error.

C. ORGANIZATION OF THE PAPER
The paper is organized as follows. The signal model for the drone RC FHSS signals is given in Section II. Mathemetical preliminaries of cyclostationarity and DF methods are laid down in Section III. Section IV presents the details of proposed method utilizing preliminaries provided in the two previous sections. The measurement setup is explained in Section V. In Section VI, measurement results are given followed by discussions. In Section VII, concluding remarks are provided and further studies are discussed.

II. SIGNAL MODEL
A single FHSS signal can be written as [19], where s(t) denotes the complex baseband equivalent of the information bearer that has a periodic burst type transmission for t ∈ [0, T ], K stands for the total number of hops during the duration of T , f k and ϕ k represent the carrier frequency and initial phase of the kth hopping, respectively. Rectangular window function, w k (t), can be expressed as, where T h is the dwell time.
The start time of the kth hop is represented as B k . Then, we define a sequence {C k }, which is the time difference between kth and (k + 1)th hop, are given as, where B k+1 and B k are the sum of all time gaps corresponding to (k + 1)th and kth hops, respectively. There are many different topologies for the antenna array. Considering the uniform linear array (ULA) enumerated with 0, 1, . . . , M elements, complex baseband equivalent of the received passband signal for any mth antenna element can be described by, where * is the convolution operator, h m (t) is the channel impulse response between the mth element and signal source, which has frequency flat fading, x(t) denotes the desired FHSS signal and n m (t) stands for the complex additive white Gaussian noise (AWGN) in which I and Q components are i.i.d ∼ N (0, σ 2 /2). Also, i m (t) represents the possible interference signal that disturbs the desired signal received by the mth element.Also, note that we assume i(t) and x(t) are uncorrelated signals. Considering the drastic increase in the number of wireless communication devices, the impact of the interference might occur as co-channel interference. The applications that use unlicensed spectrum consist of Wi-Fi technology, ZigBee signals, and unstructured signal sources VOLUME X, 2019 such as microwave ovens or several medical devices. Hence, to define the signal model for the unlicensed spectrum region (i.e. 2.4-2.48 GHz ISM band) can comprise an interference signal.

III. MATHEMATICAL PRELIMINARIES A. CYCLOSTATIONARITY
Cyclostationarity analysis can be utilized to discover hidden periodicities within a received signal [20,21]. In this study, we adopt the second-order cyclostationarity of the received signal since it reveals the hopping rate of the FHSS signal. Cyclostationarity analysis begins by taking the autocorrelation function which is denoted as, where (·) * is the complex conjugate operator, and τ represents the time lag. If the autocorrelation function has a periodicity in t, it can be written by Fourier series expansion, where A k (τ ) represents the kth coefficient at τ time lag which is also known as a cyclic autocorrelation function (CAF), and T 0 is the fundamental period. Furthermore, the frequency domain representation of the signals can extract unique features. In this regard, the Fourier transform of the CAF can be calculated by using the cyclic Wiener relationship, where S k (f ) is known as spectral correlation function (SCF) for a fixed k value.

B. MULTIPLE SIGNAL CLASSIFICATION
The well-known MUSIC algorithm is mainly dependent on the correlation matrix of the data and the ability to extract the signal and noise eigenvectors from input correlation matrix [22]. In the AAS, different time delays occur since the impinging signals arrive to the antenna elements at different times as shown in Fig. 1. Therefore, after defining a reference element, the delay between the mth element and the reference element can be calculated as, where d represents distance between adjacent elements, θ denotes the direction of arrival of a signal impinging upon the ULA, and c stands for the speed of light. The received signal at mth element can be formed by using the delay amounts as, where y 1 (t) represents the impinge signal to the reference element, ψ = − 2π λ d sin θ is the spatial frequency, and a(ψ) denotes the steering vector for the x(t). Considering all the elements in ULA and the additive noise samples, we can write the received signal in a matrix format for p signal sources as, where A is the M × p matrix to the p steering vectors, N represents additive noise for M element. The input covariance matrix calculated as, where (·) H is the Hermitian transpose operator, R xx stands for the signal correlation matrix, σ 2 denotes the noise variance with identity matrix I M . However, in practical applications, R yy usually can not be directly obtained and only sample covariance,R yy , can be used [23], where L represent the number of snapshot. In order to achieve the frequency content of the signal, MUSIC uses the eigenvalue decomposition. If the corresponding eigenvalues of theR yy are sorted in decreasing order, the largest p eigenvalues indicate the signal subspace and the remaining (M−p) eigenvalues of theR yy represents the noise subspace. Therefore, using the signal and noise subspaces, pseudospectrum of phase can be calculated as [22], where A x (θ) denotes the steering vector for x(t) and U n refers to the matrix of eigenvectors associated with the noise subspace. A search over (13) is performed to find maximum points as,θ whereθ is the estimated AoA of the signal.

C. ROOT-MUSIC
The root-MUSIC method, in comparison to the MUSIC algorithm, finds the roots of a polynomial rather than plotting the pseudo spectrum or searching for peaks in the pseudo spectrum [24]. Same as the MUSIC algorithm, the covariance matrix can be estimated from several snapshots. Then, eigenvalue decomposition of the estimated covariance matrix is employed to obtain spectral function. The polynomial can be obtained by taking the inverse of the MUSIC spectrum and is given as [25], where steering vector and noise subspace are the same as in (13). Considering the roots of the polynomial that are inside the unit circle, the closest p roots to unit circle are selected. Finally, the estimation of the AoA can be determined by, where d is the distance between two adjacent elements, λ denotes the wavelength and ψ stands for the roots of f r (θ) −1 .

IV. PROPOSED METHOD
Although FHSS signals occupy a wide spectrum, each hop of FHSS signals has a lower bandwidth. On the other hand, an accurate estimation of the AoA needs continuity in the spectrum. For the case of FHSS signals, presence in a certain frequency with a short time is not enough for an accurate estimation. Many hops are combined sufficiently, the accumulated signal samples can be increased in a certain frequency which is the pre-processing step of our proposed method. In Fig. 2 the flow graph explains how the system works in brief. First, the received multi-channel signal is achieved using the ULA and cyclostationarity analysis of the received signal obtained from the first antenna is performed. After ensuring the cyclostationarity feature detection (CFD) of the drone RC, parameters of the signal are estimated to reconstruct the FHSS signal for the received signal array.
Here, it is crucial to reconstruct the correct hops of the FHSS signal. Therefore, if there is an interference signal, outliers are removed from the parameters table according to the spectral localization of the signals. Furthermore, the resampling process is employed to decrease the signal sample rate and reduce the computational complexity. In the last step, AoA of the FHSS signal is achieved by using subspace-based algorithms.

A. SIGNAL DETECTION: CYCLOSTATIONARITY ANALYSIS
In order to decide whether a drone RC signal is present or not, CFD is adopted [20,21]. Many FHSS signals follow a pattern in time and this leads cyclic behavior in time domain. The received signal can be reformed by combining the x(t) into the r(t) by, (17) The periodicity is provided by s(t) whereas x(t) and r(t) are not periodic due to the frequency hopping pattern. To analyze hidden periodicity inside r(t), frequency hopping pattern is suppressed on r(t) which is defined as, where r(t) represents the received signal from any of the antenna elements. Considering the (17) , r(t) is expressed as, Autocorrelation function is considered for cyclostationarity analysis as (21) in which R r, r (t, τ ) represents autocorrelation of r(t). R r, r (t, τ ) has 81 summation terms and 64 of them are 0 due to properties of noise. σ 2 stands for the noise variance, and P s and P i are the average power of s(t) and i(t), respectively. In addition, if there is no interference signal then (21) is simplified as, Recall that s(t) is periodic that depends on hopping duration, which indicates that R s,s (t, τ ) is also periodic. Thus, R r, r (t, τ ) can be expanded via Fourier Series coefficients. Finally, SCF for FHSS signals can be obtained by calculating the Fourier Transform of the CAF which defined in (7). SCF should have peaks at cyclic frequencies for which the fundamental frequency is defined as a hopping rate.
Considering the ISM spectrum band, FHSS signals such as Bluetooth and drone RC signals can share the frequency channels without frequency planning. Identification of FHSS signals from other signals even though they use the same technologies is crucial to proceed with the DF. In this regard, CFD is applied to estimate the hopping rate of the FHSS signal. Since the cyclostationarity based methods generally VOLUME X, 2019 focus on periodicities such as chip rate, pilot signals, etc., non-periodic signals such as DSB-SC AM do not indicate cyclic features. Hence, any signal that behaves periodically is detected by CFD. In order to show that, among the various communication systems that use FHSS, drone RC and Bluetooth signals are chosen for the test purposes. In the context of the 2.4 GHz ISM band, Bluetooth signals are one of the most expected signals, which also uses the FHSS technology. Please note that cyclic frequency is approximately 100-200 Hz for drone RC signals according to hopping duration [26,27]. However, cyclic frequency is 1600 Hz for Bluetooth signals which uses the BLE 4.2 technology. Moreover, if the ISM spectrum is occupied with various signal sources, the CFD analysis estimates the periodicity feature for signals that are periodic and reveal peaks for them. However, the searching area of cyclic frequencies for the Futaba T8J drone RC signal is limited to 147 Hz as will be discussed in Section V [28]. Hence, the existence of other peak values would not be a concern to consider. Therefore, we first apply the CFD to identify the unknown FHSS signals and then proceed to DF estimation of the drone RC signal. Please note that the above cyclic frequencies can be off for other values.

B. RECONSTRUCTION OF FHSS SIGNALS
A signal received in a wideband spectrum bears much higher noise level for each hop of an frequency hopping (FH) signal that makes it almost impossible to distinguish it from the desired signal. Therefore, extraction each hop from the wideband spectrum and estimation of parameters of the FH signal are required. This can be achieved by performing a two-step post processing approach. First, we apply the time-frequency analysis by using STFT to estimate the time durations, center frequency, and bandwidth information for each hop. In order to evaluate the STFT, time dependent window is employed to the captured signal and discrete Fourier transform (DFT) is applied to the resulting window [29,30]. Second, we form a set which includes several temporal parameters (e.g., dwell time, start/stop time) and bandwidth of each hop. If there is a signal that is not satisfying to the statistic of the set, this hop is excluded because of a possible interference signal. A signal is omitted if its parameters are not in range of 2-σ from the mean of the set. Since the ±σ covers ≈ 68% of the normal distribution, outliers can be separated. The hops of FHSS signals must start after the previous-hop ends and also, the hopping duration and bandwidth of the hops are the same as each other. Thus, in the case of multiple drone RC signals randomly occupying the ISM band, the hops of the interested FHSS signal can be clustered. According to the aforementioned information, parameters that do not fit the 2-σ are labeled as outliers. Then, the received signal can be filtered by using the lower and upper frequency band, i.e, [f L , f H ], of each hop for a certain start and end time interval. A bandpass filter is designed to obtain each hop from the spectrum and suppress other signals in the same time interval. This approach also eliminates the unwanted noise from the signal. After filtering the hops of the FH signal, each hop is shifted to achieve equivalent baseband signal considering estimated center frequencies. In addition, the effects from the resolution of the time-frequency analysis lead to the slicing of the signal with margins. Since the eigenstructure-based algorithms separate the signal and noise into subspaces, it is expected that result of the direction of arrival is unnoticeable due to resolution errors. Letẑ(t) be the filtered and shifted version of each hop to its baseband frequency. It can be described as, where r(t) is the received signal,f k = (f L + f H )/2 denotes the estimated center frequency of kth hop, B k = k−1 l=1 C l and C l refers to the time difference between start time of lth and (l + 1)th hop.
By assuming that the sequence of time gaps is periodic with N , C l can then described in terms of the estimated time gaps as, Then,ẑ k (t), which is baseband equivalent and also shifted to the interval [0, T h ] of the kth hop, can be expressed as, In order to increase the signal resolution at baseband frequency, time gaps between hops should be removed. Therefore, we define a concatenated function,f (t), as, As a result, by combining (23), (25), and (26) the whole reconstruction process from r(t) tof (t) can be formed as: where t k ∈ [(k − 1)T h , kT h ]. Note this process is realized for every received signal from each antenna.

C. RESAMPLING
In the final step before executing the AoA algorithm, resampling process is employed. Sincef (t) is still a wideband signal, the signal can be modeled as non-frequency hopping, which eliminates the spreading and has the bandwidth of each hop. Hence, it can be converted into a relatively narrowband signal. In this regard, resampling is applied by the factor of fractional rate which is calculated according to the bandwidth of each hop of the FHSS signal and sampling rate [31,32]. Algorithm 1 demonstrates the process for an input-output relationship as a pseudo-code. The input signal is upsampled by adding zeros between samples of the original signal. After that, an FIR anti-aliasing filter was applied to eliminate discontinuities. In the last step, filtered signal samples are discarded to decimate the signal, and samples are kept at each downsample step size. It is important to note that, phase of the input signal can change while applying the resampling process [31]. However, since the resampling process shifts the phase of each reconstructed signal in the same manner, the phase difference between channels does not change. Thus, the result of the AoA is not affected. Furthermore, this approach also reduce the computational complexity during covariance matrix evaluation.

Algorithm 1: Resampling process
Input: Complex baseband signal (f (t)), estimated bandwidth (b w), sampling rate (f s ) Result: Resampled signal (f RS (t)) Initilizations: N = length(f (t)) n ← the numerator for the fractional rate betweenb w, f s d ← the denominator for the fractional rate betweenb w, f s The AoA estimation of FH signals is performed for signals captured by using over-the-air received signals by using the ULA. The measurements are taken in the test field of TUBITAK BILGEM in Gebze, Turkey. Measurements are conducted at a suburban area with a hilly terrain structure and foliage, which is close to the Sea of Marmara for both LOS and NLOS conditions, as represented in Fig. 3(a).

A. HARDWARE SETUP
The test-bed consists of the Futaba T8J drone RC as an FHSS signal source, four identical quasi Yagi antennas, and one National Instruments (NI) PXIe receiver to support the multiple input structure as shown in Fig. 3(b). In the testbed (a) Futaba T8J RC is used as an FH signal source. The signal source operates in 2.4 GHz ISM spectrum band, (b) Four identical antennas are utilized to construct our AAS. In AAS, it is preferred to utilize ULA for the AoA process. A ULA structure is constructed with four identical quasi Yagi antennas. The separation distance between each adjacent antenna element is kept as λ/2 ≈ 6.2 cm where λ is the wavelength at 2.42 GHz. Furthermore, the height of antennas VOLUME X, 2019 is set at approximately 1.5-meter and thus, reflections from the ground are avoided, (c) The signal from each antenna are received synchronously with the help of a NI PXIe-1065 four-channel receiver. NI PXIe-1065 receiver chassis consists of four major parts: one NI PXIe-8108 embedded controller, one NI PXI-5652 RF signal generator, four NI PXIe-5622 digitizer which has a 16-bit resolution and four NI PXIe-5601 RF downconverter in which each downconverter covers the frequency range 10 MHz to 6.6 GHz and has a 50 MHz instantaneous bandwidth. The corresponding block diagram of our setup is shown in Fig. 3(c).

1) Description of the Futaba T8J RC FH Signal
Futaba T8J RC signal source employs the first 50 MHz of a 2.4 GHz ISM frequency band which is divided into 30 RF channels each with a channel width of 1.5 MHz [28]. The illustration of the hopping pattern of the FH signal source is shown in Fig. 4. ∆t 1 , ∆t 2 , and ∆t 3 represent the time gaps between the hops and their repetitions. These parameters are given by, The fundamental period of the Futaba T8J RC signal can be calculated as [28], 3T h + ∆t 1 + ∆t 2 + ∆t 3 C1+C2+C3 = 6.8 ms (29) where T h defines the dwell time and ∆t 1 < ∆t 2 < ∆t 3 .

B. EXPERIMENTAL PROCEDURES
Since the signal source operates in the 2.4 GHz -2.45 GHz spectrum, it spans over 50 MHz bandwidth. At this point, while the signal is captured at a high sample rate, only a section of the spectrum has been considered to prevent data overflow. Therefore, the bandwidth of interest is adjusted to 10 MHz and resulting sampling rate is 20 MS/s. The center frequency is set to 2.42 GHz to monitor the spectrum where the hops of the FH signal will likely appear the most. The signals are collected as an I/Q samples with a length of 3 MS. Also, MATLAB 2018B is used to run the pre-processing algorithm for the FH signal, the MUSIC, and root-MUSIC algorithms.
One should note that the accuracy of the MUSIC and root-MUSIC algorithms highly depends on the phase of the received signal. If there is a phase difference between each receiver before the measurement campaign, the performance of the algorithms degrade tremendously. In order to ensure the phase coherency between each receiver, we perform a calibration procedure to cancel out the possible phase mismatches that can occur from local oscillator errors and environmental factors such as cable length etc. [33].

1) Calibration Process
First, a signal generator is configured to generate a narrowband signal at 2.42 GHz. The signal generator is placed right across the ULA to guarantee 0 • for AoA. Four-channel receiver is set at 2.42 GHz. Phase differences of each antenna are calculated according to the reference antenna. Also, this process is implemented in the LabVIEW program. Finally, once the phase difference are calculated, each phase difference value is integrated into the measurement process before starting the measurement campaign. Thus, phase coherency is assured across the receivers.

VI. MEASUREMENT RESULTS
In this study, FH signals have been reconstructed with respect to parameters of each hop, and the performance results of two different AoA algorithms for FH signals are evaluated under LOS and NLOS cases.
First, the drone RC signal detection is accomplished by applying CFD to provide distinction from other signals. In this regard, the received signal is simulated for the Futaba T8J drone RC signal with respect to measurement setup. The Monte Carlo simulation is analyzed under the Rayleigh fading model for each channel and SNR ranges between -18 dB to 2 dB. The probability of detection is obtained with cyclostationarity based signal detection for different channels of the receiver as shown in Fig. 6. Since the channel is modeled as Rayleigh fading, the detection probability varies slightly between each channel of the receiver. However, 0.9 detection accuracy is reached after the -4 dB SNR for all channels of the receiver. Hence, cyclostationarity is employed for the real-world considerations. For instance, the drone RC signal and the Bluetooth signal are identified with the maximum peak of the cyclostationarity function that described the hopping rate of FH signals as seen in Fig. 5(a) and Fig. 5(b). Since the CFD is based on matches between adjacent signals, more hops that captured will give a higher peak value in the hopping rate. Furthermore, the number of hops that had been received for each distance, which is showed as drone RC locations is 20, 24, 22, 11, and 4 for the observation time interval, respectively.
After being assured of the drone RC signal, in order to obtain reliable results a minimum number of hops should be captured by the receiver. Fig. 7(a) shows the hopping pattern of a captured signal and how the received signal behaves during the observation time. The parameters of the FH signal are extracted and reconstruction process is applied as discussed in Section IV. Fig. 7(b) shows how the sparse hopping signals are assembled together and shifted to baseband frequency to get the continuity. This is employed for the received signal from each antenna.
The DF results obtained for the received signals with preprocessing and without pre-processing stages for MUSIC and root-MUSIC algorithms are shown in Table 1, Table 2,  Table 3, Table 4, and   The drones can be moved from one point to another point within a short time. Besides the fast movement of the drone, the drone controller does not move that fast. Therefore, we demonstrate the DF of the drone RC with slow angle changes for different distances. The results are given in Table 1, Table 2, Table 3, Table 4, and Table 5 based on different antenna numbers, processing effects, and different estimators. The results show that angle accuracy is improved with antenna number and pre-processing effect. In the literature, root-MUSIC known to be manifest a convergent behavior under the assumption of increased number of antennas. Therefore, it is expected to compare MUSIC and root-MUSIC performances in a generalized measurement campaign. However, due to the number of antennas available in this campaign no significant difference is observed between MUSIC and root-MUSIC algorithm. The comparison of estimated AoAs for MUSIC and root-MUSIC algorithms illustrates that the pseudo spectrum is shifted with a miscalculated value. For instance, the phase-spectrum of the measurement that is taken from a 115.24-meter distance is calculated for with processing and without processing as shown in Fig. 8. It clearly indicates that pre-processing of the received signal array improves the DF estimation of FHSS-type drone RC. However, the estimator results are getting worse with fewer antennas (e.g. M=2) for both with pre-processing and without pre-processing.
In order to validate the proposed method, the results are compared with [17]. In [17], Anderson-Darling test goodness-of-fit approach which depends on statistical characteristics of the signal is utilized for signal detection. However, the performance of statistical methods for spectrum sensing algorithms decreases in cases where the statistical VOLUME X, 2019   characteristics of noise and fading components are dominant [34]. Also, this approach can determine the presence of the signal, but can not determine whether the detected signal belongs to the drone RC signal. On the other hand, cyclostationarity signal analysis is much robust under the fading and low SNR values [35]. Since cyclostationarity analysis takes advantage of higher-order statistics, the cyclic features could still be extracted under difficult channel conditions and variable noise levels [36,37]. Furthermore, in the proposed method drone RC signal can be classified by using CFD which results of the SCF has a cyclic frequencies at a hopping rate of the FHSS signal. Thus, in the real world scenarios CFD method performs well.
It is also noted that environment characteristics of the measurement campaign affects the performance of the proposed algorithm. In Table 5, the results are calculated for the signal that is captured at the farthest distance while considering LOS condition. According to the Table 1, Table 2, Table 3,  and Table 4 which are the NLOS cases, the increase in the Tx-Rx separation degrades the results. On the other hand, at 512m which is the LOS case, the difference between true value and estimated value is only 1.10 • . Considering sensor fusion for the drone detection suite which is provided by combination of image processing, acoustic, RADAR systems, the environment conditions can be deceptive. However, the drone's pilot signal can be detected and the direction of the signal can be estimated by passive RF sensing for both LOS and NLOS cases. This additional information will narrow the search space for other sensors such as camera of the drone detection suite.
At this point, one might wonder the performance comparison for various direction-finding algorithms that has different methodologies. In order to calculate the root-meansquare error, Monte Carlo simulation is considered for MU-SIC, pseudo-Doppler, correlative interferometer, and expectation maximization under the AWGN condition. As seen in Fig. 9, the subspace-based algorithms have prominent results among other methods. While the SNR is high the correlative interferometer and expectation-maximization methods diverge, however the error values of MUSIC and pseudo-Doppler methods decrease as SNR increases. Moreover, as discussed in Section I, the pseudo-Doppler method has insufficient resolution for direction estimation. Hence, these information leads to the usage of subspace-based algorithms.

VII. CONCLUDING REMARKS AND FUTURE DIRECTIONS
This study focuses on the problem of AoA estimation of real-world drone RC signals. For this purpose, the signals are collected over-the-air by using quasi Yagi antennas with ULA for different locations of the signal source. First, cyclostationarity based signal identification is applied to distinguish drone RC signals from Bluetooth signals. Once the drone RC signal is identified, rather than directly feeding wideband signals to the AoA algorithm, the performance of the DF with the reconstruction of FH signals is discussed. For this reason, time-frequency analysis is used to get the FH signal to reconstruct in the baseband center which leads to the accurate estimation of AoA. The performance results show that gathering the hopping signal samples at a frequency point against the noise will improve the performance of DF estimation for FH signals. In future studies, we will consider the Kalman filtering for the tracking of the hops of an FH     signals and use different subspace-based AoA estimation methods. Also, for the case that the multiple drone RC signals exist, the probability of correct clustering under the different conditions for various SNR values will be analyzed. Furthermore, if there is a high resolution DTED-2 digital map of measurement region available, footprinting based location finding can be further applied to locate the drone RC. .