Space-Time Adaptive Cancellation in Passive Radar Systems

A critical issue in the realization of passive radar systems is the effective suppression of the zero Doppler interference (ZDI). The performance of the clutter cancellation relies much on the used algorithms. Several state-of-the-art approaches consist in the independent use of spatial and temporal algorithms for ZDI suppression. In this paper, a novel interference cancellation algorithm is proposed, which jointly exploits the available information from both space and time domains. We call this novel method Space-Time Adaptive Cancellation (STAC), and it differs from previous schemes included among the tools of SpaceTime Adaptive Processing (STAP). The STAC algorithm is able to cancel out the high power direct and multipath signals on all the surveillance antenna channels prior to the beamforming and other space domain processing methods. The proposed preprocessing technique facilitates the detection of targets even in the direction of the illuminator or large static scatterers. Moreover, the interference rejection capability is also enhanced in the space domain. The performance of the presented algorithm is verified through simulations and field measurements. The experiments show improvement over the currently available filtering techniques.


Introduction
The detection performance of passive radars is strongly limited by the direct and zero Doppler multipath signals from the used illuminator of opportunity (IO).The ambiguity function of the typically available reference signals (e.g., Digital Video Broadcasting-Terrestrial (DVB-T)) often has limited dynamic range.Therefore, a target can only be detected when its corresponding correlation peak can rise above the side lobes of the reference signal components.This is the so-called masking effect that appears in the case when unwanted echoes exist in the received signal in addition to the useful target reflection [1].In recent times, the ZDI suppression is an actively researched field of passive radars.
Recently published papers present the applicability and effectiveness of spatial filtering algorithms.One of the most commonly used adaptive beamformers is the minimum variance distortionless response (MVDR) algorithm.It performs the filtering based on the knowledge of the target direction of arrival (DOA).The beampattern is evolved by these kinds of algorithms in a way that it places nulls in the direction of high power interferences that are uncorrelated with the desired signal.The operation of these algorithms is widely inspected in passive radar scenarios [2][3][4][5][6][7][8][9].
Besides the application of spatial filtering by means of digital beamforming, also the adaptive temporal cancellation methods have a significant role in ZDI suppression.These techniques in passive radar scenario utilize the available reference signal to cancel out the zero Doppler-shifted and time-delayed replicas of the transmitted illumination signal.These procedures are quite different in implementation but solve the Wiener filtering problem without exception.The resulting FIR filter produces the sum of the properly weighted and delayed replicas of the reference signal, which is then subtracted from the surveillance channel [9][10][11][12][13][14][15][16].
As the proper interference suppression has major importance in passive radar systems, the spatial and temporal filtering must be utilized simultaneously.The traditional procedure is to combine the antenna channels with the surveillance beamformer and in the second stage perform the temporal filtering.This implementation is the straightforward way.In [9], the system performs the optimal MVDR beamforming in the first stage and then filters the residual direct path interference with the CLEAN algorithm.The same filtering schema can be seen in [6,15,17,18].
The application of the filtering methods in this way may encounter some limitation as the beamforming in the first stage wastes system resources.The reference-correlated direct and multipath signals can be effectively suppressed in the time domain where the resources are quite cheap.In contrast, for the elimination of a signal component in the space domain, a full coherent antenna receiver channel is required (including antenna element, RF front-end, digitalization, and beam space processing stages, resp.).Moreover, the spatial degrees of freedom can be utilized much more favorably to the suppression of nonreference-correlated interferences.Another problem is that beamforming techniques for direct path suppression blind the passive radar systems in the direction of the IO [3,12].Avoiding blindness in the direction of high power interferences could be especially important where high azimuth coverage is required [17].Villano et al. [3] proposed the application of principal eigenvalue beamformer to avoid blind zones.This method places nulls only in the direction of the most dominant clutter components.For this reason, the suppression performance degrades.It would be rather advantageous to filter the antenna channels in the time domain first; then, the adaptive beamformer would receive a cleaner input resulting in a better distribution of the system resources and more enhanced detection.A processing scheme that performs time domain filtering on the surveillance antenna channels and then applies beamforming in the range-Doppler domain is used by the SMARP system in [5].
The implementation of this preprocessing scheme is not straightforward as the direct and independent filtering of the antenna channels (as implemented in [5]) may worsen the situation due to the fact that spatial information is disregarded.Thus, research has been carried out in this field.
The proposed algorithm is a combination of spatial and temporal filtering, using the information available from both domains to effectively cancel out the direct and multipath signals.The most outstanding feature of the proposed filtering procedure is that targets can be detected also in the direction of the high power IO.As the direct path interferences are totally cancelled prior to the beamspace processing, the applied beamformer will not place nulls in the direction of the IOs; hence, the target reflection will not be suppressed unintentionally.
Let us definitely distinguish the proposed prefiltering algorithm from the widely used STAP as the latter is mainly applied on moving platforms and operates in the angular and Doppler domain.The STAP uses the so-called spatiotemporal steering vector to match the desired signal component to the received multichannel signal matrix.In contrast, we use the term cancellation to indicate that the removed clutter components are reconstructed and subtracted at the receiver side.In Section 2, the considered space-time signal model is introduced.In Section 3, the currently used clutter cancellation techniques and the proposed prefiltering method are described.Its configuration and fast implementation techniques are also presented in detail.Finally, the presented method is validated through simulations (Section 5) and field trials (Section 6).

Space-Time Interference Model
A proper signal model is constructed at the first place.The signal propagation paths of a typical passive detection scenario are illustrated in Figure 1.The passive radar receives the emitted illuminator signal on a number of different signal paths.Besides the useful target reflection, direct and multipath signals also appear in the received signal on the surveillance antenna.In the case of Digital Video Broadcasting-Terrestrial (DVB-T) illumination, it is common to use the nearby transmitter stations in single-frequency network mode, where all the transmitter towers transmit the same signal on the same frequency with finely tuned delays.This effect complicates further the situation [7].
Let us assume that the sampled illuminator signal transmitted from one of the towers is s n , where n denotes the discrete time index.The k sample time delayed signal path has the α k complex coefficient, which describes its propagation factor.In addition, we can assign DOA information to every incident signal component.For this extension, an antenna array is assumed on the receiving side.We must also distinguish the signal components in the received signal, which has the same time delay (k) but comes from a different DOA. a ϑ k,q denotes the array response vector which belongs to the signal path that has k time delay and comes from the qth incident angle.Correspondingly, we need to introduce different propagation factors to the impinging signals that have different DOAs but have the same time delay.Thus, we use α k,q instead of α k .Figure 2 illustrates the introduced signal path notation in detail.The dashed lines show the constant bistatic range ellipsoids, while the solid lines represent the signal propagation paths.On the kth bistatic ellipsoid Q k , a number of signal paths with different DOAs are considered.In Figure 2, Q 0 = 3 and Q 1 = 2, respectively.
For an M element equidistant linear antenna system, the array response vector has the following form: where ψ k,q = 2π/λ d cos ϑ k,p , λ is the wavelength of the center frequency of the source signal, d denotes the antenna element spacing, and the direction of arrival for the impinging signal is denoted by ϑ k,q .The antenna element index is denoted by m, where 0 < m ≤ M − 1 holds.Thus, a ϑ k,q ∈ ℂ M×1 .Correspondingly, the mth channel of the received surveillance signal is described by International Journal of Antennas and Propagation In (2), f k,q denotes the Doppler frequency of the current signal component, while f s is the sampling frequency.For ZDI signal components, f k,q = 0, and for the desired target reflection component, f j,p ≠ 0 is assumed.K is the number of considered time delays, which is absolutely determined by the environment.We also consider D interference sources that are statistically independent and uncorrelated with the illuminator signal s n .Thus, E s n i d n = E s n E i d n , where E denotes the calculation of the expected value.ψ d = 2π/λ d cos θ d describes the phase shift of the dth interference component that is arriving to the array from the θ d angle.α d denotes the propagation factor of the dth interference component.The nth sample of signalindependent noise collected by the mth receiver channel is denoted by ξ m n .A more compact form of the received signal array x s n ∈ ℂ 1×M can be formulated using the array response vector a ϑ ; thus, The multichannel noise vector at the nth time instant is denoted by ξ n ∈ ℂ M×1 .In the following investigations, we assume that the illuminator signal s n and the interference signals i d n , d = 0, … , D − 1 are wide-sense stationary (WSS) processes with zero mean and unity variance.
Besides this, we can safely consider the case when the received reference signal contains only the transmitted illuminator signal plus uncorrelated noise.
x r n = s n + μ n 5 In real situations, it can be fairly guaranteed with the use of high gain antennas.Also for IOs, where digital modulations are applied such as DAB, DVB-T, DVB-S, UMTS, and LTE reference, signal reconstruction techniques are able to provide the clear reference signal [7,19].According to the described signal model, we can construct the corresponding adaptive filter.

Traditional Clutter Cancellation Algorithms
Let us first define the used reference and surveillance signal vector notations, which will be used in the outlined theoretical description of the filters.(iv) x s n is a row vector that stores the signal samples of the surveillance antenna channels at the time instant n, x s n ∈ ℂ 1×M .
(v) x r n is the reference signal vector, which contains the last recorded K samples starting from the nth signal sample, n = 0, … , N − 1.The signal samples that outlie the coherent processing interval are assumed to equal to zero.
The used notations for linear algebra operators are the following: (i) a * is the complex conjugate of the a vector.
(ii) A H is the Hermitian transpose or the conjugate transpose of the A matrix.
(iii) A ⊙ B is the Hadamard product of the two matrices.
(iv) A ⊗ B is the Kronecker product of the two matrices.
According to the presented signal model, the main objective of the clutter cancellation is to suppress all the signal components in x s n expecting the target reflection α j,p s n − j e j2π f j,p /f s n a θ j,p .
3.1.Adaptive Spatial Filtering.Adaptive beamformers combine the surveillance antenna channels using a coefficient vector in a way that the interference components are canceled and the desired signal is enhanced.The nth sample of the beam space processed surveillance signal is calculated as follows: where w s m is the mth coordinate of the beamformer coefficient vector.w s is chosen in a way to minimize clutter power.The space adaptive processing (SAP) is often realized by the optimum MVDR beamformer.In this case, w s is calculated using the spatial correlation matrix R s and the array response vector of the signal of interest a θ j,p .
Other adaptive beamforming techniques may utilize information about the eigenstructure of the spatial autocorrelation matrix [2,3].Clearly, the degrees of freedom (DOFs) of the beamformer is equal with the number of antenna elements M; hence, only M − 1 interference component can be dealt with.Also, it is important to note that no information is used regarding the illuminator signal s n .

Adaptive Temporal
Filtering.Time domain filters subtract the weighted reference signal and its time delayed replicas from the surveillance signal.In systems where the number of receiver channels is limited only, time domain filter can be used M = 1, x 0 s n → x s n .
where w t j is the jth coordinate of the temporal coefficient vector and J is the DOFs of the temporal adaptive filter.w t is calculated using the following expression [20]: where R t is the temporal autocorrelation matrix and r t is the temporal autocorrelation vector.
The calculation of the optimal w t is an active field of researches.According to the sample matrix inversion (SMI) technique, the temporal autocorrelation matrix and the cross-correlation vector are estimated with their sample average.Other techniques such as the variants of the extensive cancellation algorithm (ECA, ECA-B, and ECA-S) apply different time intervals to the coefficient estimation and filtering [10,14].Iterative algorithms like the least mean square, normalized least mean square, recursive least squares, and so on update the w t vector from sample to sample [11,16,20].Time domain filtering is a very effective method for the suppression of the time-delayed replicas of the reference signal; however, it is not able to deal with interferences that are nonreference correlated (i d n , d = 0, … , D − 1) or lie outside the time delay range of the filter s n − k , k > J.

Temporal after Spatial
Filtering.In multichannel systems, the combination of spatial and temporal domain filtering is used.The filtering process is preformed in two stages where the first stage beamformer suppresses the unwanted interference components in the space domain and in the second stage the temporal filter cancels the reference signalcorrelated ZDI components.It is illustrated in Figure 3.
This filtering procedure is used in many systems [5,15,17,18].The filtered surveillance signal for this filter scheme can be calculated as follows: where the beamformer coefficient w s * m is calculated according to Section 3.1 and the time domain filter coefficients are calculated using the algorithm described in Section 3.2 with replacing x s n to y sd n .

Space-Time Prefiltering
The main point is to realize that utilizing the spatial filtering to suppress the time-delayed replicas of the illuminator signal may result in a waste of resources.Adaptive beamforming techniques such as the MVDR try to place nulls in the direction of the interferences.It can be proved without difficulty that the beamformer response in the direction of an interference is dependent on the power of the interferer.

4
International Journal of Antennas and Propagation where a θ i is the array response vector of the particular interference, Q is the spatial correlation matrix when the interference is not present, and σ 2 i denotes the averaged power of the interference component.In most cases, the direct and multipath zero Doppler signals received from the used illuminator have very high power relative to the thermal noise floor and the target reflection Thus, the first stage spatial filter consumes the hardware resources to cancel out the high power ZDIs.At the same time, when the target direction is close to the incident angle of a clutter component θ j,p ≃ θ k,q or θ j,p ≃ θ d , the target reflection will be suppressed by the null formed into the direction of the high power interference.This behaviour is also proved with simulations presented in Section 5.
However, these interference components can be eliminated with using purely signal processing resources in the time domain.The idea behind the following presented filtering technique is to apply time domain filtering prior to the adaptive beamformer to save resources.

Spatial Extension of the Transversal
Filter.Before we introduce the optimal filter regarding the previously expressed signal model, we show a straightforward spacetime filtering scheme, which is the simplest solution but unfortunately inadequate when adaptive beamformers are used after.One can say that the surveillance antenna channels can be prefiltered by M individual temporal filter.In this implementation, the filtered output of the antenna channels is described by (15).
We can rewrite this formula into a matrix equation where W is the filter coefficient matrix in which each column w m is assigned to an antenna channel.Thus, w k,m is the weight coefficient for the k time delay replica on the mth antenna channel.The maximum assumed time delay is K.
The optimal vector for the mth surveillance antenna channel is calculated using the following formula: where r n is the cross-correlation vector of the mth channel.As it is obvious from the formula, it does not use spatial information (a θ ) in any sense to produce its output.Moreover, there is not any information coupling between the channels; thus, the signals are passed through the filter completely independently.In the optimal scenario, where there are not any undesirable RF effects on the antenna channels, this filtering scheme can work properly; however, in a real environment, it is likely to fail.This behaviour is supported by field measurement results in Section 6.3.Moscardini et al. [6] have made efforts to investigate a similar filter structure.They implemented the spatial adaptive processor in the range-Doppler domain where the individual range-Doppler matrices are obtained from the time domain filtered surveillance antenna channels.However, in their work, the spatial correlation matrix for the beamformer is estimated in a different way from the range-Doppler map.The same signal processing chain is addressed in [7].

Space-Time Adaptive Cancellation (STAC).
Following the previously described signal model (3), one can construct the corresponding filter structure to cancel the ZDI components.It is illustrated in Figure 4.The filter first prepares the timedelayed replicas of the reference signal; then each replica is extended to the space domain by applying the associated steering vectors.Finally, the proper complex weight coefficient of the current signal path is adjusted, which sets the amplitude and the common phase shift.The temporal dimension of the filter is K; thus, the clutter is expected on the first K range cells.Besides this, in each of the K ranges, Accordingly, the total number of considered interference components is The output of the filter is described by In (18), the weight coefficient of the filter is denoted by w k,q .The output of this filter can also be represented in a matrix form using Space adaptive proc.
X s (1) [n] where ψ m is a vector that contains the proper phase shift values corresponding to the direction of arrival on the mth channel and x rQ n is the extension of the previously described x r n reference signal vector where i Q k is a vector, which has Q k unity elements.
It must be emphasized that in (19) the w weight coefficient vector is common to all channels and w, x rQ n , ψ m ∈ ℂ Q×1 .In this case, all the antenna channels use the same w vector and the individual channels are linked together through the ψ m vectors.In order to perform the filtering, the optimal coefficient vector for this model must be determined.Let us define the performance function J n as the exsspected value of the squared power of the filter error.Our main goal is to find the optimal w coefficient vector for this filter architecture that minimizes the performance function.We can follow the well-known steps of the derivation of the Wiener filter to obtain the final equation for the weight vector calculation [21].where e n = x f n is defined as the error vector.In this notation, e n is a row vector.
The error vector and the vector of the filtered channel outputs are then calculated with the following matrix formula: where the matrix S Q is the extension of the reference signal vector and the Ψ matrix contains the phase compensations for the associated reference signal time-delayed replicas.This phase compensation term is related to the phased array surveillance antenna. 1) [n] X s (1) [n] X s (0) [n] Let us note that the lth row of the Ψ matrix is the array response vector corresponding to the lth signal component Ψ l, * = a θ l T .As Ψ contains information about the direction of arrival of the clutter components, we can call it as the clutter DOA matrix.
Using the above described notations, J n can be written as follows: Then, w is found by solving (28) [21].
We can write the solution in a matrix form as follows: where R Q is known, respectively, as the autocorrelation matrix and r Q as the cross-correlation vector.These terms are defined as follows: Coefficient Calculation Method.The autocorrelation matrix and the cross-correlation vector are often estimated using their sample average.

31
It is often called the sample matrix inversion (SMI) technique.The calculation of these terms is a very computationally intensive task; thus, in the following, a fast calculation technique is introduced.
Let us assume that the number of considered clutter components is a constant value for all ranges Q k = c, k = 0, … , K − 1.Note that this assumption does not restrict the operation of the filter and the further theoretical investigations as the weight coefficient corresponding to a nonexisting clutter component will go to zero when N goes to infinity.It can be proved without difficulty that R Q can be decomposed to the Hadamard product of a temporal autocorrelation matrix and a clutter DOA correlation matrix.
In this form R t Q is a block matrix, where the elements in each block correspond to the elements of the traditional autocorrelation matrix used in many time domain least squares filter.In a consequence, R t Q can be obtained very fast by simply calculating the standard temporal autocorrelation matrix and extend it to a block matrix.
where I c ∈ ℝ c×c is an all-ones matrix.
Besides this, we can develop a closed form to calculate the R Ψ Q clutter DOA autocorrelation matrix in the knowledge of the clutter DOA incident angles The ith and jth elements can be calculated as follows: where the ith row selects the clutter component with k1, q1 range and direction indexes, while the jth column selects the clutter component with k2, q2 indexes.
A fast formula can also be developed for the calculation of the r Q cross-correlation vector.Using the same principle, we can define the cross-correlation vector of the mth channel: that is, the block extension of the standard temporal crosscorrelation vector of the mth channel.
where i c ∈ ℝ c×1 is an all-ones vector.Finally, the r Q vector can be calculated as the sum of the phase compensated crosscorrelation vectors of the individual antenna channels.

39
Using the presented method, we are able to reduce the computation costs by calculating only a K × K sized autocorrelation matrix and M piece of K dimensional cross-correlation vector instead of a Q × Q and M piece of Q × 1 dimensional terms, respectively (we assumed Q = K c for simplicity).
7 International Journal of Antennas and Propagation Also note that however R , the fast formula for A ⊙ B −1 is not known; thus, the inversion of the Q × Q dimensional R Q matrix can not be avoided.
4.4.Clutter DOA Estimation.One essential question remains in the implementation of the STAC filter that is the clutter DOA information.However, one can calculate the filter coefficient using (29); the Ψ clutter DOA matrix must be a known preliminary.In practical cases, this information is not an available preliminary; thus, we must obtain it from the received signal.One way to gain knowledge of the θ k,q , k = 0, … , K − 1, q = 0, … , Q k − 1 angles is to perform several DOA estimations.It is vital to the operation of the STAC algorithm that the incident angles of the clutter components arriving from different ranges must be distinguished.It is essentially needed to prepare the S Q and Ψ matrices.We can manage to extract this information by performing DOA estimations for each range bin.In an ideal case, the surveillance signal array consists only of the following terms for the DOA estimation at the lth range.
According to the considered signal model in (3), the other unwanted signal components should not present method to approximate u l is to apply coherent integration as the s n − k , k = 0, … , K − 1 is known from the reference signal x r n .

41
Without utilizing coherent integration at this stage, the estimated DOA angles may suffer from the effects of the low SNR condition and the cochannel interferences that result in inaccurate estimation.This behaviour of the DOA estimators is deeply investigated in the literature [22].For the estimation of the DOA of the clutter components at the lth range, we can define the spatial correlation matrix as follows: It is often estimated by its sample average, After the spatial correlation matrix, R l s ∈ ℂ M×M has been obtained; the DOA of the desired signal components can be estimated using traditional DOA estimation methods like multiple signal classification (MUSIC), Capon's estimation or Burg's Maximum Entropy Method (MEM), and so on.Nevertheless, problem may arise from the fact that zero Doppler signals with the same time delay are totally correlated.
The MEM shows higher resistance against the correlated sources than Capon's or the MUSIC method.It may correctly detect the number of correlated sources, but the results will suffer from inaccuracy.Despite this fact, it is feasible in all cases and could be a reasonable choice.The power angular density with the MEM can be calculated using [23] where s θ denotes the steering vector and r j is the jth column of the inverse of the spatial correlation matrix.The MUSIC algorithm is often used in phased array systems due to its high angle resolution capability.It utilizes the eigenstructure of the spatial correlation matrix to find the direction of uncorrelated sources.The spatial spectrum can be estimated with the MUSIC using with E being the noise subspace matrix, which is spanned by the eigenvectors corresponding to the smallest eigenvalues of the R l s matrix [24].Since the basic algorithm (Schmidt-1977, Bienvenu-1979) fails to resolve correlated sources, a number of research have been carried out to resolve this weakness.Widespread preprocessing techniques for resolving the coherent source localization problem are the forward-backward averaging and the spatial smoothing [22,23,25].The spatial correlation matrix lost from its rank and became nondiagonal and singular in a coherent environment.These methods are intended to restore the full rank of the spatial correlation matrix.They rely on the fact that the array response vector a θ in uniform linear arrays is the same when its elements are reversed and conjugated.The forward-backward averaging is able to decorrelate two sources by replacing the standard spatial correlation matrix with its forward-backward smoothed version [24].This matrix is given by where J ∈ ℂ M×M is the exchange matrix, which is defined as follows: The spatial smoothing technique can be applied to resolve even more coherent sources in the expense of 8 International Journal of Antennas and Propagation decreasing aperture.Using the forward-backward spatial smoothing technique, it is possible to detect up to 2M/3 coherent sources.To achieve this improvement, the antenna array is divided into several subarrays.Let us denote the subarray size with L. Then from an M element antenna array, P = M − L + 1 subarray can be formed.The subarray formation is illustrated in Figure 5.
The signal vector of the pth forward subarray is made by the selection of the corresponding channels in the surveillance signal array.In this particular case, we use the coherent integrated u l n channel vector instead of the raw surveillance channel vector x s n .Then, the forward subarray is denoted by v f p n and can be mathematically described following formula: where p ∈ ℕ, p = 0, … , P − 1.In contrast to the forward subarray selection, we can define the pth backward subarray as Note that in the backward subarray the elements are also complexly conjugated.Now the forward-backward spatially smoothed correlation matrix can be defined as It must be emphasized that using this technique the effective aperture of the antenna system is decreased to L from M antenna elements, R l s−ss ∈ ℂ L×L [23,24].Qing and Ruolun made efforts to determine the sufficient and necessary conditions to decorrelate the signal sources with spatialsmoothing.Their work can be very useful to understand the limitations of this technique [25].
After all, the proposed overall clutter DOA estimation procedure is illustrated in Figure 6.As it was described previously at the first step, the signal of interest is enhanced with coherent integration.After that, forward-backward averaging or spatial-smoothing is preformed in order to decorrelate the reference signal components that have the same time delay.After the spatial-correlation matrix is obtained, the DOA estimation method can be used to determine the power angular density information.
In the next step, the dominant peaks must be found from the DOA method output.The peak find process generates a list of peaks that may correspond to the DOA of the sought signal components.Antenna systems that have interelement spacing larger than half-wavelength suffer from the aliasing effect [3,26,27].This results in multiple peaks corresponding to the same impinging signal.Of course, it is sufficient to select only one DOA for the filtering algorithm from among the overlapped peaks.Thus, to save resources, we may eliminate the peaks lying outside the unambiguous region.The borders of this region can be determined using The last processing section in the clutter DOA estimation is the peak selection.In this stage, the available spatial DOFs must be taken into account.As it is well known, only M − 1 interference can be cancelled using an M element antenna system.Correspondingly, M − 1 direction must be selected for the filtering and the others must be discarded even when those correspond to dominant clutter components.Thus, we can declare the following: It is worthwhile to mention, when Q k = M the clutter DOA information is totally disregarded in the filter, because the filter has enough dimension to produce spatially independent weight coefficients.In this case, we return to the spatial extension of the transversal filter, which is described in Section 4.1.The determined clutter DOAs are then stored for every k time delay, and finally, the Ψ clutter DOA matrix can be prepared for the STAC filter.

Filter Summary.
A filtering scheme that utilizes the STAC technique as a prefilter is illustrated in Figure 7.The main idea behind the presented filter is to cancel the ZDI prior to any space domain processing to save system resources.As it is highlighted in the previous section, the available degrees of freedom in the space domain limits the number of cancelled interferences Q k ≤ M − 1, ∀k 0, … , K − 1 .However, in such environments where the number of interferences at a given range is lower than M − 1, the filter is able to cancel all the disturbances and the time domain filter stage is not required.
A fast procedure for the STAC filtering can be summarized as follows: (1) Estimate the clutter DOA information and prepare the clutter DOA matrix Ψ according to the processing scheme illustrated in Figure 6.(5) Calculate the Q dimensional autocorrelation matrix R Q and cross-correlation vector r Q using (32) and (39), respectively.

Simulation
The effectiveness of the proposed algorithm is investigated by comparing its performance with the available space, time, and space-time domain filters using simulated data.On the hardware side, a 6 element uniform linear array is assumed to receive the surveillance signal, while the reference signal is assumed to be received on a dedicated antenna and receiver channel.The samples of the illuminator signal are drawn from a normal distribution having zero mean and unity variance s n ∼ N 0, 1 .The simulated clutter and target parameters are summarized in Table 1.The first set of signal components describes the considered parameters of the ZDIs, while the next one corresponds to the target reflection.Accordingly, the clutter components have zero Doppler frequency f k,q = 0 Hz.The final parameter set describes the considered jammer.The indicated α k,q and α d parameters give Clutter DOA  DOA method (MEM, MUSIC, etc.) 1) [n] 1) [n] x r [n − l], l = 0, … , K − 1 P() . . . . . .X S (1) [n] X S (0) [n]  Table 1: Simulated signal parameters.The simulated signal is generated according to formula (2).Using the previously described signal composition, a number of different effects can be investigated at the same time.The incident angle of the clutter at the 2nd range cell is 40 °which is very close to the target DOA that may result in unintentional target echo suppression.At the 3rd range, cell clutter is received from multiple directions; thus, coherent signals must be resolved by the clutter DOA estimator.At the same time, 4 directions are assumed, which is the theoretical maximum of the number of resolvable coherent signals by the spatial smoothing technique (51).

Performance Metrics.
The performance is measured with the clutter attenuation (CA) previously introduced by Cardinali et al. [16] and with the estimated target signal to interference plus noise ratio (SINR).Clutter attenuation measures the averaged power of the cancelled signal.Let us denote the processed surveillance channel with y n .Then, the clutter attenuation is defined as We also calculated the estimated target SINR, which is defined by the following term: where χ k, f is the cross-correlation of the filtered surveillance channel and the time and Doppler shifted replica of the reference signal.
Besides that, the K and ℱ sets have the indices of the neighbouring range-Doppler cells leaving out the closest 2 cells as the guard cells.L denotes the number of considered cells in the calculation.That is, the clutter power around the target peak is estimated.

Processing Setups and Results
. The simulated passive radar signals x r n , x s n m , m = 0, … , M − 1 are processed with different filter setups in order to demonstrate and highlight the STAC filter relevance.
(i) In the first case, no filtering is applied.In Figure 8(a), the calculated range-Doppler matrix can be seen.
As it is obvious from the figure, the desired target reflection at k = 25 and f = 125 Hz is totally masked by the clutter components.
(ii) In the second processing scenario, both space and time domain filterings are utilized.This method is described in Section 3.3.The steering vector 11 International Journal of Antennas and Propagation of the MVDR beamformer is set to the target DOA a θ j,p = a 41 ∘ .The calculated beampattern is shown in Figure 9 with the solid line.We can observe that the spatial filter was not able to effectively deal with all interferences as the number of clutter components is higher than the DOF of the digital beamformer.The jamming source at 130 °also could not be suppressed.The beamformer also placed a deep null in the direction of the strong clutter at 40 °, and thus, it suppresses unintentionally the target reflection at 41 °.In the second stage, the time domain filter is configured to J = 6.The calculated range-Doppler map in Figure 8(b) shows that all the ZDI components are successfully cancelled and the target can be also identified.However, the jammer is not filtered, and the clutter is filtered in the space domain in the expense of the suppression of the useful target echo.
(iii) In the third filter setup, the STAC is performed prior to the beamforming and the time domain filtering.
The filtering scheme is illustrated in Figure 7.The STAC-filtered channels x f n m , m = 0, … , M − 1 are calculated using the method described in Section 4.2.In Figure 10, the results of the clutter DOA estimation can be seen.We can observe that using the proposed technique in Section 4.4 the clutter components at different ranges could be measured and separated successfully.In the next stage, the beamformer output is calculated using (8) with replacing the raw surveillance antenna channels x s n m to its filtered versions x f n m .The calculated beampattern for the filtered surveillance channels is shown in Figure 9 with the squared line.As it can be seen, no nulls are placed in the direction of referencecorrelated clutter sources as they are cancelled prior to the beamformer.This results in proper main beam formation in the direction of the target.Besides this, the freed DOFs could be used to filter the uncorrelated jamming source at 130 °.The output of the beamformer is filtered in the time domain at the last stage.The time domain-filtered signal is obtained with using (10) and setting J = 6.In Figure 8(c), the calculated range-Doppler map can be seen.We can conclude from the figure that all the clutter components are suppressed successfully.
As it is obvious from the simulation results, the traditional filter approach fails in several ways in the considered scenario.However, these problems can be resolved by the application of the STAC method.Table 2 summarizes the obtained simulation results by applied filtering techniques.The improvements can be clearly seen from the table both in the sense of the target SINR and the clutter attenuation.The table shows results for the cases also where only the spatial filtering or the temporal filtering has been applied according to Sections 3.1 and 3.2.

Field Trials
Besides simulations, the signal processing gain of the STAC algorithm is verified with field measurements.The proposed algorithm assumes statistically independent illuminator signal samples.For this reason, samples drawn from uncorrelated normal distribution are used through the simulations.In a real environment, the DVB-T signal is very close to uniform power spectral density and thus suitable for demonstration purposes in this case.
6.1.Demonstrator Hardware.The receiver hardware used for the demonstration is a quad channel UHF band receiver especially designed for passive radar application.The detailed description of the used demonstrator system structure and its capabilities can be found in [26].One receiver channel from the four is dedicated to receive the reference signal.The remaining three channels are connected to a quad element linear antenna array with interelement spacing of d = 0 528λ, where λ is the center wavelength of the illumination signal.Leaving the third channel unconnected in the surveillance antenna system, it can be used in minimal redundancy alignment (MRA).This arrangement is applied to have maximum aperture size for DOA estimation.The antenna channel connection and the receiver build up are illustrated in Figure 11.
6.2.Scenario and Evaluation.The field measurement was taken near to the Ferenc Liszt International Airport at Budapest, where the receiver system has recorded several landing airplanes.The utilized illuminators were three DVB-T transmitters operating in single-frequency network mode at 634 MHz with horizontal polarization.The recorded samples were processed offline.In the first step, the landing airplane is detected and its range-Doppler track is extracted.Using the target track, the SINR is calculated with the different filtering algorithms to gain knowledge on their behaviour.The estimated SINR of the target is then used as the performance metric.The STAC algorithm is implemented in the way as it is illustrated in Figure 6 except the MRA modifications.After the spatial correlation, the matrix 12 International Journal of Antennas and Propagation has found that its deficient elements (R i,3 , R 3,i , i = 1, … , 4) are completed utilizing the correlation matrix Hermitian property R H s = R s .

Measurements in Real
Environment.It must be strictly emphasized, with using three antenna channels, only two coherent sources can be resolved and thus Q k ≤ 2, ∀k ∈ 0, … , K − 1 in this particular measurement.The hardware resources of the used demonstrator system are extremely limited, but the interference cancellation improvement over the currently available techniques is still evincible.
Figure 12 shows the extracted SINRs of the captured landing airplanes.The different curves correspond to different domain processing techniques.In Figure 12(a), each of the SINR curves is relative to the processing case when no filtering is applied.The gain of the different algorithms is summarized in Table 3.The investigated filtering setups are the same as described in Section 5.2.In addition, the   13 International Journal of Antennas and Propagation filtering setup that is described in Section 4.1 is also listed.In this processing setup, the standard temporal filter precedes the beamformer.Thus, the surveillance antenna channels x s n m , m = 0, … , M − 1 are filtered individually according to (15); then, the filtered output is calculated with beamforming of the prefiltered antenna channels using (8), where x s n m is replaced by x f n m .The temporal dimension is uniformly J = K = 128 in all cases, and the MVDR beamformer uses a 83 ∘ as the steering vector to receive the target reflection.
In Figures 12(b) and 12(c), the results achieved with the two different prefiltering algorithms can be seen.The SINR curves are relative to the standard space-time filter when no prefiltering is used (Section 3.3).The processing configurations are summarized in Table 4. Configurations 2, 3, and 4 belong to the STAC prefilter.These curves also illustrate the effect of the different configurations of the clutter DOA estimator.We can conclude from the measured curves that both the alias filtering and the spatial decorrelation have improvement on the clutter DOA estimation and thus on the STAC algorithm.
As it can be seen from the figures, the method described in Section 4.1 shows no improvement (configuration 1).This is due to the fact that spatial coherence is violated and the adaptive beamformer can no longer operate properly.The obtained beampatterns for a snapshot with the different prefiltering techniques are visible in Figure 13.It can be clearly seen that using the STAC prefilter the MVDR beamformer could realize a deep null to cancel the interference at 140 °; thus, it successfully contributed to the appropriate utilization of the system resources.
6.4.Results and Discussion.The field trial results presented in this paper are obtained from records of landing airplanes, in order to ensure a well-controlled measurement environment.This configuration can grant the small variation of the grazing angle of the approaching targets and also a fairly consistent SINR intesnity increase.These considerations can emphasize the efficiency of the SINR-based performance improvement evaluation.At the same time, the receiver system was only able to record and transfer the samples with roughly 60 percent of duty cycle due to the large amount of generated data.As a result, the analysed radially approaching targets are observed only in a limited number of snapshots (less than 50).This could yield a poor estimation and hence a larger uncertainty on the expected value of target SINR improvement.However, as it is visible on the obtained results, the SINR curves belong to the differently processed surveillance records which follow the same tendency (see Figure 12(a)) and the improvements appear as an offset.This can lead us to the conclusion that the improvement is clearly visible and can be applied to a trustworthy validation of the proposed algorithm.
The realized processing gain over the currently available algorithms with this system is roughly 1 dB, while the overall gain from the clutter cancellation is 14 5 dB maximum.The STAC algorithm is mainly intended to make the space domain processing methods more effective by performing early spatially coherent clutter mitigation.The achievable improvement relies much on the manageable clutter components that are directly proportional to the DOF of the surveillance antenna system.It can be concluded that the proposed algorithm could realize improvement in a real environment even with strongly limited hardware resources.In this particular field measurement, the demonstrator system was only capable to cancel two clutter components at any given range cell due to the extremely low antenna element count.According to the obtained results, the clutter had even more components that means only a fraction of the disturbances could be cancelled prior to the beamformer.This leads to a reduced overall improvement.In a more complex system, the achievable improvement increases with the extension of the surveillance antenna array with more elements until all the clutter components are totally cancelled in the space-time domain.The better distribution of the system resources that is the main contribution of the proposed algorithm could also be observed on the measurement results (see Figure 13).

Further Research Directions
The proposed filter is derived for idealized conditions.We preliminary assumed a stationary environment and statistically independent signal samples.In real cases, the illuminator signal samples may be partially correlated.Thus, to improve the filter performance, we must deal with these problems in clutter DOA estimation and filter coefficient calculation.
Besides this, the nonstationary behaviour of the environment could also decrease the performance.Using shorter time intervals than the coherent processing interval for the estimation of the clutter DOA and the filter coefficient could deal with the problem.Moreover, the clutter DOA estimation could be quite difficult, inaccurate, and computationally intensive to calculate for all coherent processing intervals; thus, it would be beneficial to learn more about its statistical behaviour.
Extension of the filter in the Doppler dimension [14] is not studied in this paper.Further works may investigate its performance.
(i) x r n is the nth signal sample of the reference signal.(ii) x m s n is the nth signal sample of the mth surveillance antenna channel.(iii) x m f n is the nth signal sample of the mth filtered surveillance antenna channel.

6Figure 1 :
Figure 1: Signal model of the passive radar scenario.
arg min w J n = arg min w e n e n H , 22

Figure 4 :
Figure 4: Structure of the proposed interference canceller filter block.

Figure 6 :
Figure 6: Signal processing chain for clutter DOA estimation.
and Propagation the amplitude of the propagation factor of the current signal component.The values are expressed in decibels relative to the direct path signal amplitude (k = 0).The phases of the propagation factors α k,q and α d are generated randomly with uniform distribution for all components.The received signal is corrupted with white noise having −50 dB power relative to the direct path signal E ξ m n 2 / α 0,0 2 = −50 d B. The coherent processing interval in the simulation is N = 2 20 sample.

Figure 8 :
Figure 8: Calculated range-Doppler matrices with different interference cancellation schemes.The dynamic range is uniformly 60 dB.

Figure 10 :
Figure 10: Clutter DOA estimation results for the simulated signal.The dynamic range is 60 dB.(a) Clutter DOA estimation on the first 4 range cells.(b) 2D map of the clutter DOA for the first 6 range cells.

Figure 11 :
Figure 11: Block diagram of the demonstrator hardware and its antenna channel connections.
STAC before space and time domain filtering

Figure 12 :
Figure 12: SINR improvements obtained from real measurements.(a) Improvements from the different domains (target number 1).(b) STAC improvement with different clutter DOA estimation techniques (target number 1).(c) STAC improvement with different clutter DOA estimation techniques (target number 2).