Separation and localization of multiple distributed wideband chirps using the fractional Fourier transform

In this paper, we consider the problem of localizing the multiple distributed wideband chirp sources using the fractional Fourier transform. The model in the time domain and that in the fractional Fourier domain derived by the Taylor series expansion are presented respectively. The representation of location vector in the Dechirping domain is illustrated which is only related to the central angle. A novel direction of arrival (DOA) estimation algorithm in the Dechirping domain is proposed, which is extended from the conventional multiple signal classification (MUSIC) algorithm in the time domain. Using this algorithm, except for estimating the DOA, the incidence source number can be determined as well which is allowed to exceed the sensor number in the array. To demonstrate the performance of proposed algorithm, numerical results are conducted. Compared to the previous FrFT-MUSIC algorithm based on the assumption of point source model, the proposed algorithm performs a better estimation performance, especially for large angular spread and low signal-to-noise ratio.


Introduction
The direction of arrival (DOA) estimation for wideband chirp signal has been greatly achieved under the point source assumption [1][2][3]. However, the signal is always affected by multipath and scattering propagation in the practical application, which will cause the angular spread, such as the local scattering source in the mobile multipath environment, the distributed target reflection wave in the low-elevation radar target tracing system, the surface and the bottom reflection signal in the sonar detection of the shallow sea, tropospheric or ionospheric propagation of radio waves, the part of detection target in the passive radar and sonar system [4][5][6]. In this situation, the received signal in an array can be considered as a superposition of scattered signals originating from the different direction and the estimation performance of DOA will be degraded significantly if using the traditional estimators based on the point source model [7,8].
For obtaining the exact nominal DOA and angular spread of a spatially distributed source, the problem of distributed source model has been widely studied since the early 1990s, and a large number of methods are proposed for the parameter estimation of distributed source [9][10][11][12][13][14][15][16][17]. However, most of the models and the estimation algorithms can be only exploited in the case of narrowband source, because the location vector in the time domain is time varying when the incident source is wideband. To date, the estimation method for wideband source with local scattering is still in scarcity. M. Ghogho et al. proposed the distributed wideband source model [18]. Liu et al. [19] proposed a kind of maximum-likelihood (ML) estimates for finite-bandwidth distributed sources by the perturbation method. Foroozan and Asif [20] introduced time reversalbased range and DOA estimators to exploit spatial/multipath diversity existing in strong multipath environments, which has a high complexity in the multi-source scenarios. Mecklenbrauker et al. [21] extended the Bayesian approach to a distributed wideband source for the application to seismic recordings, the kind that should be sparse.
The wideband chirp signal, whose frequency is linearly increasing with time, is widely used in many applications, such as communication, radar, sonar, and biomedicine, and also be used as signal model for a good deal of natural phenomena. However, for the distributed wideband chirps, we have not found much work referring to the problem of DOA estimation.
In our previous work [22], the properties of the fractional Fourier transform (FrFT), which is a generalized form of the ordinary Fourier transform (FT), is exploited to address the problem mentioned above because it provides a compact representation for chirp signals as it is based on the decomposition of the signal on the orthonormal basis set of the chirp functions [23,24]. The models of wideband chirp sources in both time domain and energy-concentrated domain are studied. And an estimation of the spatial parameters of multiple wideband chirps with local scattering is provided.
In order to obtain better estimation precision, we will further study the problem in another fractional Fourier domain, the Dechirping domain, in this paper. The representation of location vector in the Dechirping domain is derived based on a number of good properties of chirp signal in two special fractional Fourier domains, which is only related to the angular parameter. Then, the standard MUSIC algorithm combined with the source separation technique is extended to estimate the incident angles of multiple sources in the Dechirping domain. Numerical results illustrate that the proposed algorithm could separate multiple chirp signals successfully, which is far more than the number of sensors, and estimate the DOA of each chirp accurately.

Notation and definition
The fractional Fourier transform as a linear integral transform with kernel K α (u, t) [25,26]: where F α [⋅] denotes the FrFT operator and kernel K α (u, t) is where the domain u is generally known as the fractional Fourier domain which makes the rotation angle α with the time domain.
Rather than defining the fraction of the transform as the rotation angle, α, in the interval [−π, π] radians, a new variable, p, is defined as the order of the transform and is valid in the interval [−2, 2]. The zeroth order transform is simply the function itself, whereas the first order transform is its Fourier transform.

Time-frequency rotating properties
According to the time-frequency rotating properties of the FrFT, the Wigner-Ville distribution (WVD) of the FrFT of a signal can be interpreted as the coordinate rotating form of the WVD of this signal [27]. The energy spectrum of a finite chirp signal shows a fin-shape line in the fractional Fourier domain, as shown in Fig. 1 and two special fractional Fourier domains are particularly noteworthy: the Dechirping domain, whose axis u coincides with the line, and its perpendicular domain, energy-concentrated domain, in which the energy distribution of chirp signal shows an obvious peak.
We assume that a chirp signal is modeled as where β ¼ a 0 e jφ 0 is a constant, a 0 symbolizes the amplitude of the chirp signal, φ 0 is the initial phase, f 0 is the initial frequency, and μ 0 is the chirp rate. Then, the rotation angles of the chirp signal in the Dechirping domain and the energy-concentrated domain can be written as The relationship between these two rotation angles can be represented by where F α [⋅] denotes the FT operator. It can be seen that the FrFT of a chirp signal in the Dechirping domain equals to the FT of the FrFT of this signal in the energy-concentrated domain. As a result, when the chirp rate μ 0 is unknown which is satisfied in most instances, the rotation angle a d can be acquired using (6) after the determination of the rotation angle a e in the energy-concentrated domain by searching the coordination of spectral peak in twodimensional plane (α,u).

Time delay representation in Dechirping domain
According to (3), we can obtain the y(t − τ) as follows [25]: where τ is the time delay. Compared with (3), it can be seen that the chirp rate of the signal is invariable after the time delay performs on this signal. Therefore, the rotation angle of the timedelayed signal is the same as that of the original signal in the Dechirping domain.
In Dechirping domain, the FrFT of a chirp signal given by (3) can be written as: where According to the time shifting property of the FrFT, the FrFT of y(t − τ) in the same Dechirping domain is: where τ is the time delay.
In the digital calculation process in FrFT [28], the chirp rate will be significantly reduced because of the dimensional normalization. As a result, the numerical value of the rotation angle α d is far from π/2. If τ < < f 0 , which is satisfied in most practical applications, we can obtain that Thus, (11) can be approximated to Comparing (8) and (12), we can obtain: In (15), A(τ)is only related with the time delay τ , and the influence of the independent variable u has been eliminated. This useful conclusion will play an important role in deducing data model in Dechirping domain in the next section.

Data model in time domain
A scenario with Q uncorrelated sources transmitted wideband chirp signals is considered. Due to multipath propagation, each source can be seen as a superposition of N q scattered point-source components. The uniform linear array (ULA) with P sensors is taken as an example for derivation, and the similar conclusions are easy to be generalized to the arrays with other types. The complex envelope of the output vector in the array can be modeled as where β qi : the complex amplitude of the ith scattered source from the qth source; a θ i þθ qi ; t : the time-varying location vector of the array, θ q is the nominal DOA of the qth source and θ i þθ qi is the DOA of the ith scattered components from the qth source; s q (t): the qth transmitted signal; v qi : the time delay associated with the ith scattered source from the qth source; n(t): the additive zero-mean noise vector, which is assumed to be spatially white and independent of the transmitted signals.
When the received signals are wideband and with local scattering, there is no one-to-one relationship between the DOA parameter and the location vector a θ i þθ qi ; t which is also dependent on the time and deviation angle variables. Therefore, traditional parameter estimation methods cannot work well in the case of distributed wideband chirp. The time-frequency properties of chirp signal in the two fractional Fourier domain mentioned above can help us eliminate the influence of time and deviation angle variable on the location vector. The detailed derivation will be given in the next section.

Approximating the spatial signature
Because the majority of scattered signals distribute in the vicinity of the transmitted source, the time delay v qi is relatively small so that (12) is satisfied and the quadratic term of v qi can be ignored. Hence, the ith scattered component of the qth source received on the conference sensor in the Dechirping domain is given by Then, the FrFT with a rotation angle α d performs on the output vector in (16), and the result can be given by where A q θ q þθ qi is the location vector of the qth source, S α d q u ð Þ is the FrFT of sq (t) and N(u) is the FrFT of n(t) in the Dechirping domain.
Because FrFT of multiple sources or that of one single source's multi-scattering components does not result in cross-term interference, the conclusion in [31] could be exploited here directly as follows: where T is the transpose and where τ p ¼ p−1 ð Þd sin θ q þθ qi =c represents the time delay on the pth sensor, d is the inter-sensor distance of ULA, and c is the transmission speed. In order to estimate θ q , [29] proposed an approximate model called the generalized array manifold (GAM) and a corresponding algorithm using a Vandermonde structure. Due to the deviation angleθ qi of each scattered signal is relatively small, A q θ q þθ qi in (20) can be approximated through the first-order Taylor series expansion: Then, we define a variable v q , that is where c iθ qi and we assume that Hence, the compact matrix notation V can be given by Finally, the model (18) can be approximated as (24), we know that the spatial parameters θ À and ϕ À are needed to estimate in the approximated model.

Single source estimation algorithm
Firstly, we consider one chirp source (Q = 1 in (16)). Many DOA estimation algorithms for the narrowband point source model can also be applied here because the location matrix in the Dechirping domain has a one-to-one relationship with the spatial parameters. In this paper, the classical MUSIC algorithm in the time domain is extended to estimate spatial parameter of wideband chirp source with local scattering. The spatial estimation spectrum of standard MUSIC algorithm is In the Dechirping domain, the calibrated location vector is exploited to replace the location vector in (28), that is where Ā(θ) = [A(θ) ϕd(θ)] and ϕ θ ð Þ ¼ 1 ϕ ½ T . Hence, the estimation spectrum for the proposed model should be revised to The parameters θ and ϕ can be determined by the spectral peak position.

Separation of different chirp sources and the source number determination
As (5) shows, chirp signals with different chirp rates have different rotation angles in their energyconcentrated domain in which chirp signal shows an obvious spectrum peak. The energy-concentrated domain of each signal can be determined by the two dimension search in the time-frequency plane (α,u). In each energy-concentrated domain, only one chirp signal can acquire the best energy-concentrated property while the energy distributions of the other signals and noise are dispersing, therefore, this chirp signal can be separated from others through a band-pass filter. Then, the parameter estimation algorithm proposed above can be applied to filtered single source after that is transformed into the corresponding Dechirping domain. The detailed separation process could be found in [31]. It is noted that the separation process cannot only resist the interference of noise and other signals, but also solve the problem of source number restriction in many traditional DOA estimators. For example, the MUSIC algorithm cannot be applied when the source number exceeds the number of the sensors in the array. However, the proposed algorithm avoids this problem through the separation processing. The parameters of each source can be estimated respectively.

Simulation study and results
To demonstrate the performance of the proposed algorithm, an ULA with eight sensors is employed to estimate the spatial parameters of the proposed model. In the simulation, two uncorrelated wideband chirp sources are considered as signal emitters. The initial parameters of the first wideband chirp source are f 1 = 100 Hz and μ 1 = 10 Hz/s. And that of the second one are f 2 = 200 Hz and μ 2 = −10 Hz/s. Besides, the central angle of the first source is taken as θ 1 = 30 ∘ with the extension width Δ 1 = 2 ∘ from a uniform angular distribution, and the central angle of the second source is taken as θ 2 = − 45 ∘ with the extension width Δ 2 = 2 ∘ from a Gaussian angular distribution of width, respectively. The signal-to-noise ratio (SNR) is set at 10 dB, and the number of snapshots is 500.
The two-dimensional searching results of two sources on the reference sensor are shown in Fig. 2. According to the coordinates of two distinct spectrum peaks, the rotation angles of their energy-concentrated domain can be determined, and the source number can be easily obtained by the peak number. In Fig. 3a, b, the energy spectrums of two sources on the first sensor in the corresponding energy-concentrated domain show an obvious energy concentration respectively. In this energy-concentrated domain, the chirp signal forms a sinc function and the majority of the energy spectrum focuses on its support, and the single peak will appear on the other sensors as well. Therefore, the desired chirp signal on each sensor can be separated conveniently using the band-pass filtering, and the filtered results are illustrated in Fig. 3c, d.
The FrFT-MUSIC algorithm proposed in [30], which is applied to estimate the DOA of the wideband chirp source based on the point source model, and the method in [22] is chosen as the comparing algorithm to demonstrate the performance of the proposed algorithm for the wideband chirp source with local scattering. For each figure, 100 Monte a b Carlo simulations were run to estimate the root meansquare error (RMSE) of the estimates. We focus on the estimation of the DOA, while the other parameters are considered mainly as nuisance parameters. Figure 4 shows the RMSE for the DOA estimates versus the different SNRs using three algorithms. For both sources, the RMSE of the proposed algorithm is less than that for the other two algorithms. Figure 5 clearly demonstrates that the proposed estimator has excellent performance versus the different extension widths. With the increasing extension width, the proposed algorithm has better location accuracy than the others.
Finally, the resolution performance of the proposed algorithm versus the different source number is demonstrated. In this experiment, the chirp rate interval Δμ of each source is chosen as 1 Hz/s. The numerical results can be seen in Fig. 6. The proposed algorithm can remain the good resolution performance even though the source number is far more than the sensor number.

Conclusions
In this paper, parametric localization of multiple wideband chirp sources with local scattering have been considered. The models in both the time and fractional Fourier domain were proposed. In the Dechirping domain, a novel DOA estimator combined the properties of the chirp signal and the Taylor series expansion was addressed. The DOA of each source and the source number can be determined using the proposed method. Besides, the source number is allowed to exceed the sensor number in the array. The simulation study demonstrated that the proposed algorithm provided a superior spatial resolution performance for wideband chirp sources with local scattering than the previous FrFT-MUSIC algorithm.