Abstract
For the purpose of simultaneously estimating and locating the linear frequency modulation (LFM) emitter with unknown parameters, an innovative approach, i.e., Fast Direct Position Determination (FDPD), is introduced. The proposed algorithm is based on maximum likelihood estimation (MLE) and spectrum detection. To improve the accuracy and overcome the dramatic complexity of plain maximum likelihood formulation, we further derive the objective function equation of Direct Position Determination (DPD) algorithm and present an enhanced strategy to solve the highly nonlinear optimization problem. By combining the two-step localization method, one-step localization method, and short-time Fourier transform (STFT), our approach realizes jointly estimation of the transmitted signal parameters and emitter localization. Simulation results show that the proposed method is superior compared to the existing DPD, and two-step localization algorithms in terms of localization error and computational complexity, especially for low signal-to-noise ratio (SNR).
1. Introduction
Localization of emitters is a fundamental demand in many application fields, such as radar [1, 2], sonar [3, 4], global navigation satellite systems (GNSS) [5, 6], wireless sensor network [7, 8], and seismic exploration [9, 10]. Linear frequency modulation (LFM) signals, known as large time-bandwidth product and high resolution signals, are frequently used in localization systems and have attracted increasing attention of many scholars. As stealth technology and low intercept technology get prosperous and mature, passive emitter localization technology has gradually developed from the traditional single station to using several widely distributed receiver stations. This paper considers a passive localization problem for an emitter by multiple diffusely separated receiver stations.
In general, the localization methods contain two main categories of two-step localization method and one-step localization method. Algorithms of two-step localization method [11–20] are commonly composed of two steps: first, each receiver station measures and computes the parameters related to the position of the target. The parameters involved in position include time difference of arrival (TDOA) [11–13], angle of arrival (AOA) [14–17], frequency difference of arrival (FDOA) [18–20], etc. The second step is to construct equations using the parameters to estimate the position of the target. Both steps of two-step method might be optimal; given their input, two-step method is still suboptimal. Since the correlation between different receiver stations is ignored, the precision of the parameter measurement in the first step is limited. Moreover, the parameters may be virtually absent in some of the receiver stations due to the signal strength fluctuation and different path attenuation, leading to data association and performance deterioration issues, especially in low signal-to-noise ratio (SNR) scenarios.
One-step localization method [21–35] processes the observed signals in the central processor without parameter extraction and then constructs the cost function only related to the position of the emitter. This kind of localization algorithms finally obtain the position of the emitter through the maximum likelihood estimation (MLE), the least squares, the grid search, gradient-based methods, etc. Therefore, it is also referred to direct position determination (DPD) [21–24]. In 2004, Weiss first proposed direct position determination technology [21], analyzed, and compared the DPD performance with AOA for signals with known (DPD-known) and unknown waveforms (DPD-unknown) [22, 23]. Over the past dozen years, multiple DPD algorithms have been presented to enhance performance of passive localization system based on different approaches including the maximum likelihood (ML) [25, 26], the time frequency analysis (DPD-STFT-w) [27], the multiple signal classification (MUSIC) [28, 29], the expectation maximization (EM-DPD) [30, 31], the minimum variance distortionless response (MVDR) [32–34], and the time-varying delay [35]. Compared with the one-step approaches, the DPD algorithms directly estimate the initial position of the target according to the received signals in the central station, avoiding the process of correlation of parameters of different receivers. Moreover, DPD can locate the signals (simultaneous signals with same frequency, etc.) which are difficult to be processed by traditional methods. At the same time, because the cost function of DPD is only related to the position of the emitter, it has higher localization accuracy and better robustness at lower SNR on account of making full use of the prior information of the signals which are coming from the same emitter [23]. However, since the DPD algorithms deal with the original sampled signals rather than the analytic solutions of the emitter, the computational complexity increases rapidly and the estimation performance of waveform parameters deteriorates drastically when the SNR is lower than a threshold point, i.e., the “threshold effect.”
This paper addresses passive localization of an LFM emitter with unknown signal parameters. Aiming at the problem of “threshold effect” and high computational complexity of DPD algorithms at low SNR, we proposed a Fast Direct Position Determination (FDPD) algorithm. The proposed algorithm can realize the combination of emitter localization and parameter estimation of the LFM signal. The main procedures of this approach include the following: (i)Firstly, the signal spectrum is obtained through short-time Fourier transform (STFT) and other operations (linear fitting, the least square method, etc.); then, the SNR is calculated(ii)Compared the SNR with the threshold point, the cooperative localization of two-step method and one-step method is combined to finally achieve the optimal localization. The performance of the two-step method and one-step method is similar at high SNR, we choose the two-step method owing to its less computational complexity. Chan algorithm [36] and Taylor algorithm [37], which are two representative algorithms of two-step localization method, are selected to complete collaborative localization: Chan algorithm is used to estimate the initial position coordinates, regarding it as the initial value, and Taylor series expansion algorithm is used to obtain the further position coordinates(iii)At low SNR, the initial frequency, chirp rate, and other signal parameters are estimated by time frequency analysis, thus we can capture the entire parameter set of the LFM profile of the signal and obtain an objective function. Then, a four dimensional search (4D) is carried out at the position coordinates obtained by Chan and Taylor algorithms. To solve the highly nonlinear optimization problem (four-dimensional search), we turn the four-dimensional search into a two-dimensional parameter estimation combined with grid search. The search range is determined by the error analysis of Chan and Taylor algorithms. Simulation results show that the proposed algorithm reduces the computational complexity, and the positioning accuracy is better than Chan, Taylor, DPD, and other traditional algorithms
The rest of this paper is organized as follows. Section 2 introduces the model of the LFM signal for widely separated receiver stations. In Section 3, the localization problem formulated in an ML estimation paradigm is clarified. Section 4 presents the proposed FDPD algorithm in detail. Section 5 provides the numerical simulations. The conclusions are drawn in Section 6.
2. Signal Model
Consider a 2D Cartesian coordinate system where a noncooperative stationary emitter transmitting signals. We let the vectors represent the positions of L widely separated receiver stations. The transmitter position is denoted by the vector of coordinates . As shown in Figure 1, the emitter transmitted an LFM signal . Then, the signal is intercepted by the receiver stations and transmitted to the center processor.
The LFM signal is modeled as where denotes the initial frequency of the LFM signal, is the chirp rate, and is the time duration. The signal observed by the th receiver station is given by
is a time-dependent vector, is an unknown complex scalar representing the attenuation coefficient related with regard to the path from the emitter to the receiver stations, and is the th array response to signal transmitted from position . is the signal waveform, transmitted at time and delayed by . The vector represents the independent and identically distributed zero-mean complex Gaussian white noise with corvariance matrix , and each receiving channel is independent of each other. can be calculated by [35]
with denoting the speed of light.
After sampling, the continuous signal of Equation (2) can be expressed in a vector form where
with denotes the number of the samples, denotes the sampling frequency, is the signal bandwidth, and is the time duration of the received signal, thus ().
3. Problem Formulation about Direct Position Determination
Derived from [38], the MLE of the unknown parameter vector can be estimated by examining the likelihood ratio, with corresponding to the target-present hypothesis. When the distribution characteristics of noise are known, MLE can be used [21, 39]. Here, we assume that is the identically distributed zero-mean complex Gaussian white noise, with correlation matrix . In addition, the receiver stations are distributed widely enough to warrant mutual independence of noise vectors,
According to Equation (4), we can construct the probability density function (PDF) as where is a constant, and represents the operation of conjugation transpose. Thus, the likelihood ratio function can be obtained as
Let , , , and . Finally, we can get the maximum likelihood function as follows:
Remarkably, we observe from Equation (9) that the parameters to be estimated in the objective function include the transmitter position, the signal parameters, and the time . It is a highly nonlinear optimization problem whose closed-form description of the analytic optimum solution is difficult to obtain. The numerical solution based on grid search for the region of interest is an alterative. However, the high-dimensional grid search will become computationally prohibitive when the needed localization accuracy increases. If the SNR is lower than a threshold point, the estimation errors rise quickly, i.e., the “threshold effect.” Based on the above arguments, a novel algorithm is investigated in the subsequent sections to trade off algorithm performance for implementation complexity.
4. Algorithm Description
In this section, we present a computationally efficient approach to solve a highly nonlinear optimization problem in passive localization. First, we employ STFT technique to extract the valid spectrum fragment of the transmitting signals and then further estimate the instantaneous frequency, spatial power, and other signal parameters. After extracting peak points of the power spectrum, aiming at the disadvantages of low precision of two-step method and high complexity of DPD algorithm at low SNR, we adopt corresponding objective functions to extract the emitter positions. The algorithm is described in detail below.
4.1. Estimation of Signal Parameters and Spectrum Detection
Under the circumstance of low SNR, the localization performance of DPD algorithm will deteriorate significantly. In order to improve the localization performance, in this subsection, we adopt the time-frequency analysis method to estimate the relevant parameters and detect the spectrum.
4.1.1. Parameter Estimation Based on Time Frequency Analysis
To estimate the LFM signal parameters, the most important step is to estimate instantaneous frequency (IF) from the received signal. The IF expression of the LFM signal is as follows:
After sampling, can be written in a vector form as where , and is the number of samples of the LFM signal with .
As can be seen from the above equations, IF is formed by linear combination of initial frequency and chirp rate . A large number of IF samples can be obtained through estimation, and the corresponding and can be estimated from the IF samples by means of parameter fitting.
For the LFM signal, short-time Fourier transform (STFT) performs analysis on the signal in the time window and can intuitively reflect the trend of the frequency and the change of spectrum along with time, thus the initial frequency and chirp rate of the signal can be estimated efficiently, which has a good effect on the estimation accuracy, especially at low SNR [40–45]. The common STFT is defined as where is the time window function, which is called analysis window as well.
Substituting Equation (12) into the expression of the LFM signal [42], then the STFT for the LFM signal is defined as where
By taking the square of the module of , the spectrum (SP) of the LFM signal can be written in the following form:
The Gaussian window is chosen as the analysis window so that is a Gaussian function in essence. is symmetry about the line and achieves the maximum value at . The above formula shows that the LFM signal’s spectrum obtains the maximum value at frequency point at each discrete time . Therefore, the IF of the LFM signal can be acquired by finding the max point at each discrete time .
Hence, we can get the IF estimation as the following expression:
Figure 2 shows the real value and estimated value of the instantaneous frequency after linearly fitting by the least square method when . Thus, the estimated values of chirp rate and initial frequency can be calculated from the fitted straight line.
4.1.2. Extraction of the SP of the LFM Signal
The SP of the observed signal of the th receiver station can be expressed as [38] where represents the SP of noise.
After adding the moving window function , the extracted LFM signal can be written as
The SP of the LFM signal observed and extracted by the th receiver station in Equation (17) is shown in Figure 3.
On the basis of the above process, the 2D distribution of the chirp signal’s STFT is converted into 1D relationship between frequency and time, which displays a coarse line in the time-frequency plane. It can be seen from Figure 3 that when SNR is at -5 dB, the spectral intensity of the LFM signal is much greater than that of the noise. It can be clearly seen from the figure that the time-frequency characteristics of the LFM signal reflect the sum of power spectral in each window. Summing the maximum value of spectrum of each window helps to estimate the parameters of the LFM signal and realize the spectrum detection. The localization algorithms after spectrum detection will be introduced in the subsequent sections.
4.2. The Localization Algorithm Based on Spectrum Detection
According to the “threshold effect” of DPD algorithm, the localization accuracy will deteriorate sharply when the SNR is lower than -5 dB [21]. The mean (or root mean square) of miss distance will also increase rapidly. At high SNR, DPD, Chan, Taylor algorithms have similar performance in localization, but when the SNR is lower than the threshold point, DPD methods have been demonstrated to outperform two-step methods. For most of the DPD algorithms require grid search in the final, which increases computation complexity and transmission demand in the case of lacking prior information. Therefore, Chan and Taylor algorithms are first used for coarse location. Then, grid search is carried out on the basis of the results of the coarse location, which can greatly reduce the workload. The complexity analysis will be detailed in the next subsection. The whole FDPD algorithm process is as follows.
4.2.1. Localization at High SNR
According to [46], in passive multistation localization, after obtaining TDOA measurements by Chan algorithm [13, 36], the distance differences between the emitter and the base stations are obtained, thus multiple TDOA measurements form a set of estimated equations about the emitter. We assume is the estimated position of the emitter and is the known position of the th receiver station, then the emitter position can be calculated by the following formula: where
It has been shown that the iterative algorithms (Standard least squares, hyperbolic least squares, etc.) require a definition of the initial parameters, which can significantly affect solution accuracy and computational time [47]. Therefore, to ensure the solution accuracy, the position calculated by Chan algorithm is selected for initial position of MS for Taylor expansion of the Taylor algorithm [37]. Then, the weighted least square (WLS) method [48, 49] is used to obtain the least squares estimation of Z as follows: where represents the transpose operator, and is the covariance matrix of the measurement value of TDOA. The value of the TDOA measurements covariance matrix is estimated under the assumption that there are large numbers of samples, or equivalently, long observation time. Therefore, the Cramer-Rao matrix bound (CRMB) of the estimated TDOA vector is used to represent the value of the matrix [50]. where is the bandwidth of the signal and is the trace of matrix . is the signal power spectrum, is the noise power spectral matrix, is the lower right by partition of the matrix , and is a vector of unity which has the same size as .
The selection of the initial value of iterative Taylor-series expansion method has a great influence on the location result. If the initial value is not appropriate, the Taylor algorithm will not converge. On the other hand, Chan algorithm has fast operation speed and low computation complexity. Although the Chan algorithm reduces the localization accuracy in poor channel environment, the localization results can still reflect the general characteristics between the position of the radiation source and the measured value, which is conducive to the convergence of Taylor algorithm. Therefore, we choose Chan algorithm to calculate the TDOA measured values, and the localization results of Chan algorithm are taken as the initial value of Taylor series expansion algorithm, then several iterations are carried out.
Hence, this paper adopts the collaborative localization of Chan and Taylor algorithms using Equations (19)-(22) at high SNR, which can not only ensure the small error of the initial value of the Taylor expansion, so as to ensure the localization accuracy, but also greatly reduce the amount of calculation.
4.2.2. Localization at Low SNR
When the SNR is low, the DPD algorithm is used on the basis of the location results of Chan and Taylor algorithms. To get more accurate waveform parameters, STFT is used to assist signal parameter estimation. DPD algorithm can obtain better localization accuracy on the basis of more accurate signal parameters estimation. In order to obtain more accurate positions, here, an objective function is constructed for grid search.
The IF of the LFM signal is estimated from Equation (16), so as to estimate the initial frequency and chirp rate, and then substituted into the following function, which describes the fitting error per unit time at different positions, window width, and starting time, the function is given by where , is the number of samples. is referred to the minimum cost function, thus denotes the mean minimum cost function, i.e., the fitting error per unit time at different positions. From Equation (15), represents the total likelihood ratio function of the signal, which reflects the sum of spectral power in all windows. reflects the summation of the maximum SP within the sliding window. To obtain higher parameter estimation accuracy, the following objective function is defined:
Compared with the function in Equation (9), the objective function includes more factors of indicating the parameter estimation by means of STFT, which is maximized by . Then through grid search, we can estimate the emitter position , initial time , and signal pulse width more accurately. In most instances, four-dimensional search cannot be considered as a viable solution to solve a highly nonlinear optimization problem. Therefore, we used an efficient decoupled strategy to solve the problem. The basic idea of the strategy is to turn the four-dimensional search into a two-dimensional parameter estimation combined with grid search. Specifically, based on STFT, a two-dimensional estimator is firstly designed to extract the initial frequency and chirp rate of the LFM signal in the time-frequency domain. Then, a four-dimensional search, including the previously estimated waveform parameters, is carried out, and the position of the emitter, the corresponding transmitted time, time duration, and initial frequency can be obtained by the search.
For determining when and which method should be used, that is, to decide whether the emitter signal is in high or low SNR environment, Equation (15) in 4.1.1 can be utilized. By means of maximizing the spectrum at each window, it can be seen that the maximum spectrum of the LFM signal is obviously different from the maximum spectrum of the noise. The spectrum difference between the LFM signal and the noise can be calculated and then compared with the threshold point to determine whether the emitter signal is at high or low SNR. The simulation diagram of the maximum spectrum value changing over time is shown in Figure 4.
The general pseudocodes of the proposed FDPD algorithm are shown in Algorithm 1.
| ||||||||||||||||||||||||||||||||||||||||||||||||||||||||
5. Simulation Results
In this section, simulations are brought forward to verify the performance of the proposed algorithm. Four receiver stations used in the simulation are located at (0,0)m (m denotes meter), (6000,0)m, (6000,6000)m, and (2400,6000)m, respectively. The approximate coordinates of the radiation source are (6000,4500)m. We perform extensive Monte Carlo simulations. In the following analysis, the results are collected and processed by averaging M Monte Carlo experiments. Simulation parameters are shown in Table 1, where is the number of Monte Carlo experiments.
For each Monte Carlo experiment, we assume that the error between the estimated position of the emitter and the actual position is calculated as follows: where is the estimation of of the th Monte Carlo trial.
The mean error (ME) of the estimated position is obtained by the statistics of M Monte Carlo experiments. ME is defined as
The root mean square error (RMSE) of the estimated parameters is defined as where denotes the parameter to be estimated, such as chirp rate and initial frequency.
5.1. Analysis of Simulation Results
In Figure 3, we extract the 1D frequency-time curve of the LFM signal based on STFT and perform STFT on the transmitted signal, the spectrum of which is symmetry about. In order to get the SNR for spectrum detection, according to in Equation (18), the amount of maximum SP within the sliding window can be viewed as a factor of SP. The simulation of the variation of the maximum SP in each STFT window with time is shown in the following figure.
It can be seen from Figure 4 that when the SNR is -15 dB, the maximum SP value of the LFM signal is obviously different from the maximum SP value of the noise signal. Then, the SP differential is calculated and compared with the threshold point to determine the signal environment.
Figures 5 and 6 are simulation comparison diagrams of FDPD and common direct localization algorithms.
As can be seen from Figures 5 and 6, the DPD-known algorithm can acquire the optimal localization performance for the case that the transmitted signal is completely known to the receivers. Therefore, the performance of DPD-known can be viewed as an upper bound for other algorithms. Moreover, in high SNR environment, namely, higher than 0 dB, all the four algorithms achieve superior localization accuracy and perform nearly equally. When the SNR is below the threshold point, the ME and RMSE of these algorithms rise quickly as SNR decreases. As SNR decreases, the accuracy of the DPD-unknown deteriorates rapidly. When the emitter signal is unknown, FDPD and DPD-STFT-w algorithm has better localization accuracy compared with DPD-unknown algorithm, which can ensure the accuracy that the RMSE is less than 100 m when the SNR is above -5 dB. However, considering the extremely high computational complexity of DPD algorithms (DPD-STFT, EM-DPD, etc.), Chan and Taylor algorithms are used for coarse localization, and parameters are estimated on the basis of them, then the objective function is maximized by using the estimated parameters which are more accurate. This will greatly reduce the amount of calculation, so it is called fast DPD (FDPD) algorithm. The complexity analysis of FDPD algorithm and DPD, DPD-known and DPD-STFT-w algorithms will be described in detail in the next section.
5.2. Analysis of Complexity
The complexity of Chan and Taylor algorithms is , where is the characteristic dimension. In this paper, is equals to 2. denotes the number of receiver stations. When is 4, the complexity is , which is relatively low under the simulation background of this paper.
The complexity of DPD and DPD-STFT-w algorithm is composed of the complexity of two-dimensional search, complexity of maximum likelihood estimation, and least square iterative fitting. The complexity of maximum likelihood estimation least square iterative fitting is and , respectively, both of which increase with the number of receiver stations and the number of samples . Hence, the overall computational complexity of DPD algorithm is , where and denote the number of meshes on the -axis and -axis, respectively. Compared with DPD algorithm, DPD-STFT-w algorithm needs to take STFT of the observed signal from every receiver station, whose computation complexity is . Therefore, the overall computational complexity of DPD-STFT-w algorithm is , where is the window width of STFT, is the gliding step length of STFT, is the total number of samples of the LFM signal, and is the computational complexity of the estimation of .
For the complexity of two-dimensional search is reduced significantly, then the complexity of FDPD algorithm is , where , denotes the number of meshes on the -axis and -axis, respectively. It can be seen from the complexity of the three localization algorithms that the complexity is determined by two-dimensional grid search, thus the proposed FDPD method reduces the meshes on the axes through combining parameter estimation, time frequency analysis, STFT, and two-step method with DPD algorithm.
The complexity of these algorithms is summarized in Table 2.
It can be seen from Table 2 that the complexity of FDPD algorithm has been greatly reduced under the condition of similar positioning accuracy. The combination of coarse grid and fine grid is considered in our FDPD algorithm to realize the balance between the number of grid points and the localization accuracy. Therefore, compared with DPD, DPD-STFT-w, and most of the direct localization algorithms, it can save more than of the time.
6. Conclusions
In this paper, we consider the passive localization of an LFM signal emitter with unknown signal parameters. The problem is formulated as a highly nonlinear optimization problem at the beginning; then, the FDPD approach for jointly estimating the position of the emitter and its transmitted signal parameters is presented. The FDPD method combines two-step method (Chan algorithm, Taylor algorithm) and one-step method (DPD algorithm) with short-time Fourier transform to ensure the localization accuracy and greatly reduces the computational complexity compared with DPD algorithm. Spectrum detection is aimed at improving the “threshold effect” of DPD algorithm. When the SNR is higher than the threshold point, it can achieve almost the same localization performance as the DPD algorithm when the signal parameters are known. When the SNR is lower than the threshold point, the computational complexity can be saved more than . Simulation results show that the proposed method is obviously superior to the traditional localization algorithms. However, the FDPD method is currently only supported in terms of several types of signals and some specific conditions. Further work is underway to consider the boundary conditions and the scenario with multiple emitters transmitting complex signals, including NLFM, BPSK, and QPSK.
Data Availability
The data used to support the findings of this study are available from the corresponding author upon request.
Conflicts of Interest
The authors declare no conflict of interest.
Authors’ Contributions
H.R and R.J.G derived the proposed method. H.R and R.J.G conceived and designed the experiments. H.R performed the simulations. H.R analyzed the results. H.R wrote the paper. H.J.L and R.J.G reviewed the paper. All authors have read and agreed to the published version of the manuscript.
Acknowledgments
This research was funded by the National Natural Science Foundation of China (Grant number 61671304).