Research on Direction Finding Method under Impulsive Noise Based on Nonuniform Linear Array

. Direction of arrival (DOA) estimation under impulsive noise has always been an important research area in array signal processing. The traditional methods under impulsive noise mostly rely on prior parameters and have high computational complexity. Based on the ﬁ ltering theory, we present an effective pretreatment ﬁ ltering technology to cut out the impulse mixed in the array received data and employ the nonuniform linear array to improve the estimation performance further. First, according to the amplitude characteristics of impulse noise, the pretreatment ﬁ ltering technology is proposed to cut out the impulse based on the median ﬁ lter and sliding average ﬁ lter, which is valid for both strong and weak impulsive noise. Second, the minimum redundant array is adopted to carry out array virtual expansion so that the array aperture can be increased and the estimation performance can be improved. Finally, based on the idea of matrix reconstruction, we propose the improved estimation of signal parameters via rotational invariance techniques algorithm and an improved root multiple signal classi ﬁ cation algorithm for DOA estimation. Theoretical analysis and simulation results show that the proposed method has a simple processing process, small calculation load, good array expansion ability, and excellent noise adaptability. Moreover, the proposed methods greatly improve the direction-ﬁ nding performance under the condition of low signal-to-noise ratio and strong impulsive noise.


Introduction
Direction of arrival (DOA) estimation is also known as spatial spectrum estimation or direction finding, which is an important research direction of array signal processing.It has been widely used in electronic warfare [1], radar [2], communication [3], navigation [4], acoustic communication [5], and other fields and has greatly attracted the research interest of scholars at home and abroad.A series of super-resolution directionfinding algorithms with excellent performance have emerged [6] during the past 50 years.Classical direction-finding algorithms, such as the multiple signal classification (MUSIC) algorithm and estimation of signal parameters via rotational invariance techniques (ESPRIT) algorithm [7], which are based on the subspace decomposition, have the advantages of excellent estimation performance and moderate computation.However, due to the limitation of the noise environment and array structure, the classical direction-finding algorithms cannot be applied directly in practical applications.In many practical application fields, noise does not meet the ideal Gaussian distribution but shows different degrees of impulsive characteristics, such as cosmic electromagnetic wave, low-frequency atmospheric noise, and various human interference noise, which is called impulsive noise [8].Besides, the nonuniform linear array (ULA) structure will also limit the direct application of the classical direction-finding algorithms.
The symmetric alpha-stable (SαS) distribution process can be used to model the impulsive noise, and the characteristic exponent α ð0<α ≤ 2Þ: can be used to describe the intensity of noise [9,10].The smaller the value of α, the stronger the intensity of the noise.In the early stage, the DOA estimation method based on the fractional low-order statistics was mainly adopted for direction finding in impulsive noise, such as the fractional low-order moment and the fractional low-order covariance algorithm [11,12], which can solve the problem that there is no effective second moment in the background of impulsive noise.However, such methods have disadvantages of bad adaptability to strong impulsive noise, high computational complexity and dependence on prior parameters.Based on the idea of matrix reconstruction (MR), the reconstructed fractional lower order covariance (RFLOC) method [13] was proposed, which improves the DOA estimation performance under strong impulsive noise effectively.However, it is still limited by computational complexity and prior parameters.Then, Wang et al. [14] proposed a new method based on sparse representation and iterative hard threshold optimization algorithm, which is more robust and has good resolution and estimation accuracy.The signal detection in unknown noise is addressed utilizing the interference cancelation before detection [15] method, which has a small computational load and is effective with training data containing interference.In literature [16], based on the quantum-inspired gray wolf optimization algorithm, a direction-finding method for impulsive noise background is proposed, which has good estimation performance.However, the methods based on the optimization algorithms need multidimensional optimization and have a huge amount of computation.Based on the median filter theory, a filter preprocessing technology is designed in literature [17], which can effectively weaken the impulsive noise and make the second-moment-based direction-finding algorithms can be directly applied.
Most of the classical direction-finding algorithms are based on equidistant ULA, which has good estimation performance.But it strictly depends on the Van-der-monde structure characteristics of the array.However, in many practical engineering applications, such as aircraft and ships, the uniform array cannot be arranged ideally because of space limitations.Therefore, direction finding using a nonuniform array is an important research direction, and the sparse characteristics of the nonuniform array can be utilized to improve the estimation performance [18].Sparse array can expand the array aperture and improve the estimation performance [19], such as minimum redundant array [20], coprime array [21,22], etc.For some non-ULAs with special structures, the sparse characteristics and special position relationships between array elements can be utilized, and the matrix construction [23] and other methods can be used to expand the array aperture virtually.It can be seen that the estimation performance can be improved by exploiting the sparse characteristics of the array.A closed-form solution is proposed for multiparameter estimation [24] of the mmWave polarized massive multiple-input multipleoutput, which is computational efficient and suitable for arbitrary sensor arrays.
In this paper, we propose a direction-finding method for impulsive noise based on the minimum redundant sparse array and filter processing technology, which breaks through the bottleneck of direction finding under strong impulsive noise and has excellent estimation performance in both weak impulsive noise and Gaussian noise background.First, the minimal redundant linear array (MRLA) is adopted for mathematical modeling.Second, based on the analysis of the amplitude characteristics of the impulse noise, a new filter pretreatment technology of the array received data is proposed using the ideas of median filter and sliding average filter (SAF), which can weaken the impulse noise effectively and has excellent processing effect under strong impulse noise, weak impulse noise, and Gaussian noise.Thirdly, based on the idea of MR, the array covariance matrix data is used to expand the array aperture and translate the array into ULA virtually so that the classical direction-finding algorithms can be applied.Finally, the classical ESPRIT algorithm and root MUSIC algorithm are selected for DOA estimation.Simulation results show that the proposed method can extend the aperture of the array effectively, and the calculation is simple.It has excellent estimation performance and adaptability for impulse noise.
The rest of this paper is arranged as follows: In the second part, we give out the data model based on the minimum redundant linear array and introduce the α stable distribution.The third part mainly introduces the filter pretreatment technology and the improved ESPRIT algorithm, as well as the improved root MUSIC algorithm using MR.The fourth part presents simulations and an analysis of the proposed method.The fifth part summarizes the work of this paper.

Data Model
The minimum redundancy linear array is utilized for a data model.Suppose K signals with the wavelength λ are incident on the MRLA, the number of array sensors is M, and the spacing between the 1st and the mth sensor is d m .The noise is additive impulsive noise, which is modeled with the symmetric α stable distribution.Then, the receiving data of the mth array sensor at lth snapshot could be formulated as follows: where a mk ¼ expð − j2πd m cos θ k =λÞ: , is the time delay between the 1st and the mth sensor, cosð * Þ: is the cosine function, θ k is the incident angle of the kth signal s k ðtÞ: , and n m ðtÞ: represents the noise at the mth sensor.Let YðlÞ: ¼ ½y 1 ðlÞ; ⋯; y M ðlÞ T , SðlÞ: ¼ ½s 1 ðlÞ; ⋯; s K ðlÞ T and NðlÞ: ¼ ½n 1 ðlÞ; ⋯; n M ðlÞ T , where ½ * T means the transpose of a matrix or vector, then we have the following: where A ¼ ½aðθ 1 Þ; ⋯; aðθ K Þ: , and aðθ i Þ: ¼ ½a 1k ; a 2k ; ⋯; a Mk T , is the array manifold matrix.Since α stable distribution is suitable for modeling the impulsive noise, it has received extensive attention in signal processing.The characteristic function is employed to describe the α stable distribution.
where expð * Þ: represents an exponential function, α 2 ð0; 2, ε and γ represent the characteristic exponent, position 2 Journal of Sensors parameter and dispersion coefficient of the α stable distribution, respectively.The intensity of the impulsive noise is inversely proportional to α. Signal information has many characteristics, and the statistical moment of a signal is the best.The α stable distribution only exists fractional lower order moments with 0<p<α ≤ 2. Consequently, the theory and methods based on the fractional lower order statistics are doomed to be the best choice for signal processing problems with α stable distribution.The classical methods, such as the robust covariation (ROC), the fractional lower order moment (FLOM), and the fractional lower order covariance (FLOC), are formulated as follows, where r * ij represent the fractional lower order statistics between the ith and jth array sensors.
where 0<p ≤ 0:5α.We can learn from Formula (4) to ( 6) that they all require a prior parameter p.Moreover, the value of p will influence the results of the DOA estimation.Besides, the computational complexity is huge compared with the second-order moment methods.Therefore, it is necessary to research on the second-order moment DOA estimation methods for impulsive noise.

The DOA Estimation Based on Filter Pretreatment Technology and MRLA under Impulse Noise
3.1.The Filter Pretreatment Technology for Impulse Noise.It can be learned from the property of the SαS process that the probability density function has a large tail.The smaller the value of characteristic exponent α is, the larger the tail will be.
When the impulse noise satisfies the SαS distribution, it has a much higher probability of having a large amplitude than that of the Gaussian distribution, and it shows a certain degree of impulsive characteristics.The smaller the value of α, the higher the probability of the noise producing a large impact.That is to say, the intensity of the impact is inversely proportional to α.However, the impulse of the noise is random.It only generates a huge sample amplitude in some positions.While in most positions, the amplitude is relatively quite small, but we cannot predict when the huge impulse will occur.The traditional DOA estimation methods are mainly based on fractional lower-order statistics, which require prior parameters and have a large computation load.The preprocessing methods based on the median filter can remove the impulse noise effectively, and have good estimation performance at strong impulse.Inspired by the idea of preprocessing, we propose a new pretreatment technology in this paper.
The main methodology of the proposed technology is amplitude-limiting and sliding-filtering (AS).As the signals are buried in impulse, the amplitude limiting will be of great importance in cutting out most parts of the impulse noise, and a proper threshold is quite significant.ε is the symmetric point of the probability density function in SαS distribution.When α 2 ð1; 2, ε represents the mean value.When α 2 ð0; 1, ε represents the median value.According to the property of the SαS distribution, most data amplitude is close to the mean or median value, and the amplitude of impulse is much bigger than the median or mean value.Accordingly, it is reasonable to use the mean or median as the threshold.In order to adapt to the strong impulse, the median value is chosen as the amplitude limiting threshold.
Assume the number of snapshots is L, the array received data can be expressed as Y ¼ ½Yð1Þ; ⋯; YðLÞ: .The filtering pretreatment is implemented by two steps.First, perform the amplitude limiting on the array received data, and the amplitude threshold is the median value of the array received data.The array data at the mth sensor and lth snapshot is as follows: where Y ðm; lÞ: is the absolute value of Yðm; lÞ: , TðmÞ: ¼ medianðY ðm;:ÞÞ: , and the function medianð * Þ: represents finding the median of the vector, and angleð * Þ: means finding the phase of a complex number.Obviously, the majority impulse of the noise can be cut out by the amplitude limiting.Secondly, in order to smoothing the data noise, furthermore, apply SAF to Y A , and assume the length of the window is W.
The filtered array data can be expressed as follows: where f SAF ð * ; WÞ: means performing SAF with the window length W. So far, the proposed pretreatment technology has been implemented.The impulsive noise is eliminated mostly, and the array data has an effective second-order moment.
where N AS is the residual noise.The pretreatment technology does not need to compute the fractional lower order statistics or prior parameters about α.The process is simple and convenient, which is suitable for engineering application and realization.Simulation is conducted to illustrate the property in smoothing the impulse of the presented pretreatment Journal of Sensors technology.Figure 1 gives out the processing effect via α, and the generalized signal to noise (GSNR) is 5 dB.The preprocessing method in [16] is recorded as filtering preprocessing (FP), and the proposed method in this paper is called AS for short.Figure 1 illustrates that when the signals are buried in impulse noise, the strength of the impulse increases as α decreases, and the processing effect improves as α increases.The original array received data cannot be used for DOA estimation directly.Both the FP method and the AS method can cut out most of the impulse effectively.The FP method may cause a flat top, which may lose some useful signals.
Relatively, the proposed AS method has better smoothing performance while eliminating the impulse noise.Figure 1(a) shows that when the impulse is extraordinarily strong when α ¼ 0:3, the FP method can cut out the strong impulse but may lose some detailed information, while the AS method can keep more detailed information.Figure 1(b) shows the smoothing effect improves as the impulse becomes weak, while the FP method will still cause a flat top of the signal.

The MR Method
Based on the MRLA.After the AS pretreatment to the original array data, the majority of impulse is wiped out.The rest part of the noise can be treated as Gaussian noise, which means Y AS has an effective second or higher order moment.We can obtain the covariance matrix R of the MRLA from Y AS , where ½ * H represents the conjugate transpose of a matrix, . Since the MRLA is a sparse array, the steering vector matrix A is not a Van-der-monde array, and the classical ESPRIT and root-MUSIC algorithm can't be applied directly.Therefore, we should construct a Van-der-monde steering vector matrix A v virtually by MR first of all.
As shown in the previous data model, the sensor location of the MRLA is D ¼ fd 1 ; d 2 ; ⋯; d M g: , and d m is integral multiples of d, where d is a distance equivalent to a half wavelength.According to the property of the MRLA, the M elements number of the MRLA can be expanded to PðM <P ≤ ðM − 1ÞM=2Þ: virtually.Meaning we can get a virtual ULA with P sensors located at D v ¼ f0; 1d; 2d; ⋯; ðP − 1Þdg: from the MRLA.The covariance matrix of the virtual ULA is a Van-der-monde array, which can be reconstructed from R.
R is an M × M matrix, and the element at the ith row and jth column is recorded as Rði; jÞ: , sk is the power of the kth signal; when j ¼ i, σ 2 nij ¼ σ 2 , and when j ≠ i, σ 2 nij ¼ 0, σ 2 is the power of residual noise.As D v is exactly the array structure of a P elements ULA, we can reconstruct the covariance matrix of the virtual ULA R v from R according to the Van-der-monde structure of the covariance matrix.Apparently, R v is a P × P matrix, and the element at the pth row and qth column can be expressed as follows: where 1 ≤ q ≤ p ≤ P, when p ¼ q, σ 2 npq ¼ σ 2 , and when p ≠ q, σ 2 npq ¼ 0. And where A v is an P × K array, a v ðθ k Þ: is a vector with P elements, and u k ¼ expð − j2πd cos θ k =λÞ: .Obviously, A v is a Van-der-monde matrix, and The MR method is summarized as follows.We can get a vector r v with P elements from R, and the tth (0 ≤ t <P) element is as follows: where n is the number of ði; jÞ: pair satisfying d i − d j ¼ td.
Based on the Hermite property of R v ; we have R v ðq; pÞ: ¼ conjðR v ðp; qÞÞ: , and conjð * Þ: is the conjugate of a complex number.Construct R v from r v , Formula (17) utilizes the Hermite property and the Vander-monde property of the ULA, where 1 ≤ q 1 ; p 1 ; q 2 ; p 2 ≤ P.So far, we've expanded the M elements MRLA to a P elements ULA virtually, and the covariance matrix R v of the virtual ULA can be obtained by Formulas ( 16) and ( 17).

Journal of Sensors 5
The traditional direction-finding algorithms, such as the ESPRIT and the root MUSIC, can be applied.

DOA Estimation Based on the ESPRIT and Root MUSIC
Algorithm.Conduct eigen-decomposition to the covariance matrix R v of the virtual ULA, where U and Λ are composed of corresponding eigenvectors and eigenvalues, U S and U N are the eigenvectors relevant to the signal and the noise, so as to Λ S and Λ N .According to the characteristics of the subspace, we have the following: where T is a K × K full rank matrix.
The ESPRIT algorithm is summarized as follows.Formula ( 21) can be recorded as follows: where U S1 , A v1 and U S2 , A v2 are composed of the first and last P − 1 rows of U S and A v , respectively.And according to the Van-der-monde property of A v , Based on Formula (23), U S2 ¼ A v2 T, U S1 ¼ A v1 T, then we have the following: Apparently, the eigenvalues of Q contains the information of DOAs (u k ) we pursuit, which can be calculated from: where arccosð * Þ: is the arc-cosine function, and k ¼ 1; ⋯; K.The root-MUSIC algorithm is summarized as follows.Formula ( 22) can be recorded as follows: where u ¼ expð − j2πd cos θ=λÞ: .Obviously, the order of the polynomial f ðuÞ: is 2ðP − 1Þ: , and there are ðP − 1Þ: pair roots.Each pair of roots is conjugate to each other, and K of them are distributed on the unit circle.According to Formula (29), the K roots distributed on the unit circle contain the DOA information of the impinging signals, and the DOAs can be estimated from the following formula: where b u k is the kth root distributed on the unit circle.

The Cramer-Rao Bound (CRB).
The CRB of the DOA can be represented as follows: where Γ ¼ ½θ: , b Γ i is the estimation of Γ i , J Γ is the Fisher information matrix, and where lnð * Þ: is the natural logarithm function, ∂ð * Þ ∂Γ i and ∂ð * Þ ∂Γ j means taking the partial derivative with respect to Γ i and Γ j , f Y ðYÞ: is the probability distribution function (PDF) of array data.Let GðΓ; tÞ: ¼ YðtÞ: − NðtÞ: ¼ ASðtÞ: , GðΓ; tÞ: ¼ ½g 1 ðΓ; tÞ; g 2 ðΓ; tÞ; …; g M ðΓ; tÞ: .The PDF can be expressed as follows: where f ð * Þ: is the PDF of the impulsive noise.Assume the PDF is circularly symmetric, according to literature [15], then we have the following: where f ðεÞ εdε, and the CRB can be formulated as follows:

Simulation
In this part, simulations and analysis are conducted to illustrate the properties of the proposed methods compared with the existing methods.According to the previous description, the proposed two methods can be termed AS-MR-ESPRIT and AS-MR-rootMUSIC.The RFLOC methods in document [12] and the FP-MR-ESPRIT method in document [16] are selected for comparison because of their representativeness and excellent performance in strong impulsive noise.Moreover, the combination of the proposed MR technique and the RFLOC method is implemented for the DOA estimation of MRLA, termed MR-RFLOC.If not specified, an array with eight sensors is employed in the following simulations, and the sensors are located at ½01d 2d 11d 15d 18d 21d 23d: , where d ¼ 0:5λ.Obviously, the virtual ULA of the MRLA has P ¼ 24 sensors.The generalized signal to noise ratio ρ GSNR is expressed as follows: and E½jsðtÞj 2 : is the power of the signal, lgð * Þ: is the logarithmic function with base 10.The result is regarded as right when the error is smaller than 1 °.And Formula (33) is utilized to compute the root mean square error (RMSE).
Figure 2 shows that both the two methods have similar estimation properties and can estimate 10 signals correctly under impulse noise while the number of array sensors is 8.It demonstrates that the AS pretreatment method can deal with the impulse noise properly, and the proposed MR technique can expand the array aperture efficiently.Moreover, we can learn from Figure 2 that the proposed methods are effective at different intensities of impulse noise, and the accuracy of estimation improves with the decrease of impact the noise (as α increase).Figure 2(a) shows the validity of the proposed methods when α ¼ 0:6, the impulse of noise is quite strong.Figure 2(b) shows the validity of the proposed methods when α ¼ 1:0, the impulse of noise is generally strong.Figure 2(c) shows the validity of the proposed methods when α ¼ 1:6, the impulse of noise is relatively weak, which indicates the effectiveness of the proposed methods in both strong and weak impulse noise.
Simulation B: Suppose two signals are injected into the array; the DOA of one signal is 90°, and the other is ð90 þ xÞ: °.The resolution property of the proposed methods and the contrast methods is shown in Figure 3.The GSNR is 20 dB, the snapshots are 200, the characteristics exponent α ¼ 1:0, and the times of independent repeat experiment is 500.
We can learn from Figure 3 that the estimation performance improves with the increase of the angle interval, and the proposed two methods have the best resolution property.When the angle interval is smaller than 1°, both the AS-MR-ESPRIT method and the AS-MR-rootMUSIC method may not be able to estimate sometimes.When the angle interval is bigger than 1°, the proposed two methods have outstanding resolution properties.Besides, the proposed two methods have absolutely better resolution over the FLOC method and the FP-MR-ESPRIT method.Thanks to the MR technology proposed in this paper, which is effective in expanding array aperture, improving resolution property, and improving the estimation performance greatly.Since the MR-FLOC method combines the MR technique and the FLOC method together, and it has similar resolution properties with the proposed two methods, which illustrate the effectiveness of the MR technique in improving resolution further.
Simulation C: Assume two signals are injected into the MRLA from 40°and 55°; this simulation illustrates the property via GSNR of the proposed two methods, the other three methods, and the CRB.Figures 4 and 5 show the performance curves for α ¼ 0:6 and α ¼ 1:6, respectively.The snapshots are 200, and the times of independent repeat experiments are 500.
Simulation results illustrate that the estimated property improves as GSNR increases for all five methods.According to the AS pretreatment method, the impulse noise is eliminated effectively, and the proposed two methods have the best performance whether the impulse is strong or not.Furthermore, Figures 4(a prove the stability and effectiveness of the proposed methods in impulsive noise.Figure 4 shows that the FLOC method does not work well when the impulse is strong (α ¼ 0:6), meaning that the FLOC method has relatively poor adaptivity to strong impulse noise.Moreover, the combination of the MR technique and the FLOC method can improve the estimation performance most of the time, while it might fail to work sometimes, and the RMSE gets a large value, as shown in Figure 4(a).Figure 5 verifies that the estimated property improves as the intensity of impulse decreases for

8
Journal of Sensors all five methods.The AS-MR-rootMUSIC method is slightly better than the AS-MR-ESPRIT method when the GSNR is low, and the advantage becomes less obvious as the GSNR increases.Simulation D: There are two signals injected into the array from 70°and 80°; this simulation illustrates the estimation property of the five methods via α, as shown in Figure 6.The snapshots are 200, the GSNR is 10 dB, and the times of independent repeat experiment is 500.
Figure 6 illustrates the proposed AS-MR-ESPRIT method and the AS-MR-rootMUSIC method have outstanding performance at different characteristic exponent α.The probability of the right estimation of the proposed two methods is 100% when α>0:1, showing the excellent adaptability to the  Journal of Sensors intensity of the impulse noise.It verifies that the proposed AS pretreatment method is quite efficient in eliminating both strong and normal impulse noise.Besides, the proposed MR technique can utilize the array aperture of the MRLA efficiently, and the estimation property can be improved further.Although the FP-MR-ESPRIT method can eliminate the impulse noise effectively, it cannot utilize the array aperture efficiently, and the estimation property is inferior to the proposed methods, and the property of all five methods improves with the reduction of the impulse.Figure 6(a) shows the RMSE curve of the MR-FLOC method is not smooth when the impulse is relatively strong, since the FLOC method is not suitable for DOA estimation under strong impulse noise.

Conclusions
In this paper, we present the AS-MR-ESPRIT method and the AS-MR-rootMUSIC method for DOA estimation under impulse noise based on the MRLA.The proposed two methods break through the bottleneck of direction finding under strong impulsive noise and have excellent estimation performance in both impulsive noise and Gaussian noise backgrounds.Two key techniques, the AS pretreatment method and the MR technique, are first proposed.By studying the amplitude characteristics of the impulse noise, the AS pretreatment method is proposed, which can eliminate the impulse noise effectively and has excellent adaptability for strong impulse noise.The MR technique is adopted to utilize the array aperture of the MRLA efficiently so that the array can be virtually expanded to a ULA, and the estimation property can be improved further.Theoretical analysis and simulation results show that the proposed method can extend the aperture of the array effectively, and the calculation is simple.It has excellent estimation performance and adaptability to impulsive noise.Moreover, the proposed two methods greatly improve the direction-finding performance under the condition of low signal-to-noise ratio and strong impulsive noise, which are quite suitable for engineering implementation.
matrices based on vector, Reð * Þ: means getting the real part.
) and 5(a) show that the accuracy of the proposed two methods is much closer to the CRB.The probability of right estimation shown in Figures 4(b) and 5(b)

FIGURE 6 :
FIGURE 6: The estimation property via α: (a) RMSE of estimation; (b) probability of right estimation.