# An OPMA for Robust Mutual Coupling Coefficients Estimation of URA with Single Snapshot in MIMO HF Sky-Wave Radar

Authors: Yuguan Hou; Qingguo Jin; Shaochuan Wu; Zhuoming Li
Journal: International Journal of Antennas and Propagation
 (2016)

---

## Abstract

Due to the fluctuation of the signal-to-noise ratio (SNR) and the single snapshot case in the MIMO HF sky-wave radar system, the accuracy of the online estimation of the mutual coupling coefficients matrix of the uniform rectangle array (URA) might be degraded by the classical approach, especially in the case of low SNR. In this paper, an Online Particle Mean-Shift Approach (OPMA) is proposed, which is to get a relatively more effective estimation of the mutual coupling coefficients matrix with the low SNR. Firstly, the spatial smoothing technique combined with the MUSIC algorithm of URA is introduced for the DOA estimation of the multiple targets in the case of single snapshot which are taken as coherent sources. Then, based on the idea of the particle filter, the online particles with a moderate computational complexity are used to generate some different estimation results. Finally, the mean-shift algorithm is applied to get a more robust estimate of the equivalent mutual coupling coefficients matrix. The simulation results demonstrate the validity of the proposed approach in terms of the success probability, the statistics of bias, and the variance. The proposed approach is more robust and more accurate than the other two approaches.


---

## Body

## 1. Introduction

Multiple-input and multiple-output (MIMO) radar systems, characterized by multiple antenna elements at the transmitter and the receiver, have demonstrated the great potential for increased ability of the target detection [1–7]. One important application of MIMO radar system is the HF radar [2–5], including the HF sky-wave radar and the HF ground-wave radar. However, one important challenge in MIMO systems is that the mutual coupling becomes particularly significant as the antenna element spacing is decreased. In many practical problems, the direction of arrival (DOA) is the significant information, while at the same time many classical DOA estimation algorithms suffer from sensitivity to mutual coupling, such as MUSIC algorithm [8]. Hence, the DOA estimation for the multiple narrowband signals has been a classical problem in the array signal processing. High-resolution array-processing algorithms for source localization are known to be sensitive to errors in the model for the sensor-array spatial response. In particular, unknown gain, phase, and mutual coupling as well as errors in the sensor positions can seriously degrade the performances of array-processing algorithms. With few exceptions, high-resolution source localization algorithms require an exact characterization of array, including knowledge of the sensor positions, sensor gain/phase response, mutual coupling, and receiver equipment effects. All such information is inevitably subject to errors [9]. The presence of mutual coupling distorts the phase vectors of radiation sources and the eigenstructure of the covariance matrix [10]. In [11], the effects of the mutual coupling on the direction finding accuracy of a linear array with dipole elements were studied. It was found that a known coupling did not affect the estimation performance much. In the case of an unknown mutual coupling, the performance of most high-resolution direction finding algorithms could be degraded. Therefore, a relatively precise coupling coefficient could well establish the high-resolution DOA.In recent years, a great effort has been seen in the algorithms of the calibration and the coupling coefficients estimation [12–25]. The paper [12] described a calibration algorithm that estimated the calibration matrix consisting of the unknown gain, phase, and mutual coupling coefficients as well as the sensor positions by using a set of calibration sources in known locations, which was based on a maximum likelihood approach. However, this method requires a set of calibration sources at known locations. Two methods to compensate the unknown mutual coupling were proposed in [13], while, just the same as [12], calibration sources were required in both of them. Unlike previous array calibration methods, literature [14] proposed an algorithm that was able to calibrate the array parameters without the prior knowledge of the array manifold. Literature [15] presented a new array calibration procedure for over-the-horizon (OTH) radar, using disparate sources. The method in [16] used the signals impinging on the array to carry out both the DOA estimation and the array calibration simultaneously. In order to express the coupling coefficients, the array coupling matrix was investigated in [17, 18]. Literature [17] proposed a robust subspace-based DOA estimation algorithm in the presence of mutual coupling for ULA, which was based on the banded symmetric Toeplitz matrix model. An accurate estimate of mutual coupling matrix could be achieved simultaneously for array calibration. The method of moments (MoM) was used in [18] to evaluate the mutual coupling between the elements of a given array. The MoM admittance matrix was then used to eliminate the effects of mutual coupling. Reducing the effect of mutual coupling is also an important method for high-resolution DOA estimation [19–21]. The mutual coupling effect of a compact uniform circular array (UCA) is shown in [19] in order not to affect the favorable characteristics of the FFT-based preprocessing technique but only results in a modulation of the signal component at the receiver with a diagonal matrix. In [20], the author divided the antenna into two or more identical subarrays and discussed what kind of errors mutual coupling introduced to the accuracy of the ESPRIT algorithm. In [21], the author estimated the angles when mutual coupling was significant with dummy elements. It has been shown in the pioneering work in [22] that, by applying a group of auxiliary array elements, the MUSIC algorithm can be adopted directly for DOA estimation for ULA at first. Once preliminary DOAs are estimated, the coupling coefficients can be estimated. Then, given the estimated coupling information, the full antenna array with an enlarged array output vector can be processed to refine the DOA estimation. In [23], the 2D DOA estimation in the presence of mutual coupling was presented by setting the sensors on the boundary of the URA as auxiliary sensors. In [24, 25], the modeling and estimation of mutual coupling in a uniform linear array of dipoles were discussed and a method of the mutual coupling compensation using subspace fitting was presented.However, in MIMO HF sky-wave radar, the application of DOA estimation is often in the domain of the two-dimensional range-Doppler spectra, the peak of which is taken as the single snapshot for the receive array signal processing. The classical approaches mentioned above might be degraded in the case of the single snapshot and the low SNR. In this paper, an online particle mean-shift approach (OPMA) for the robust antenna mutual coupling coefficients estimation with the single snapshot is proposed. The goal is to find an optimal or a suboptimal estimation for the mutual coupling coefficients matrix in the case of the single snapshot. Due to the fact that the estimation of the unknown mutual coupling coefficients is related to the additive noise, the estimated mutual coupling coefficients matrix is modeled as a posteriori probability density function. Moreover, a tracking process based on the idea of particle filter and the mean-shift algorithm is applied to search the most rapid increase of the probability density function.This paper is organized as follows. In Section2, the problem and the system model are briefly formulated. The spatial smoothing technique of the uniform rectangle array is introduced. Moreover, the equivalent equation of the mutual coupling coefficients matrix is derived. An OPMA is utilized in Section 3 to further improve the precision of the mutual coupling coefficients matrix estimation. Simulation results for different scenarios are shown in Section 4 and conclusions will follow in Section 5.


## 2. Problem Formulation and System Model

Consider a collocated MIMO HF sky-wave radar system with the assumption of the narrowband waveform [2]. The transmit antenna is a uniform linear array. The receive antenna is a uniform rectangle array (URA) which is shown in Figure 1.Figure 1
Geometry of the receive antenna array.Let the transmit antenna array consist ofMt antenna elements. The receive antenna array is located at the xoy plane. Assume there are Mr receive antenna elements along the x-axis and Nr receive antenna elements along y-axis. The spacing of two adjacent elements in the transmit array is d. The center wavelength is λ. The spacing between adjacent rows of the receive array is dx and the spacing between adjacent columns of the receive array is dy, and we have dx=dy=d. Assume the matrix Xr(t) is the receive data matrix of the time t, the element of which is the receive data of the ((i-1)d, (j-1)d) antenna element denoted as Xijr(t). Let τk denote the group delay (range) of the kth target (k=0,2,…,K-1), where K is the number of the targets. Also, let θkd denote the direction of departure (DOD) of elevation from the transmit antenna to the target, and ϕkaa and ϕkae denote the direction of arrival (DOA) of azimuth and elevation from the target to the receive antenna. Assume C is the mutual coupling coefficients matrix of the receive antenna array and the mutual coupling of the transmit antenna array is ignored. Let x(t) denote the Nr×Mr-variate time t-dependent vector of the received data matrix Xr(t) [2], and we have(1)xt=vecXrt=CAst+ξt,where A is the receive antenna steering matrix and vec· is a vectorization operator; we have(2)xt=X11rt,…,X1Nrrt,X21rt,…,X2Nrrt,…,XMr1rt,…,XMrNrrtT,A=aϕ1aa,ϕ1ae,aϕ2aa,ϕ2ae,…,aϕKaa,ϕKae,st=y1t,y2t,…,yKtT,ξt=ξ11t,…,ξ1Nrt,ξ21t,…,ξ2Nrt,…,ξMr1t,…,ξMrNrtT,aϕaa,ϕae=axϕaa,ϕae⊗ayϕaa,ϕae,axϕaa,ϕae=1,βxϕaa,ϕae,…,βxMr-1ϕaa,ϕaeT,βxϕaa,ϕae=exp⁡j2πdλsin⁡ϕaecos⁡ϕaa,ayϕaa,ϕae=1,βyϕaa,ϕae,…,βyNr-1ϕaa,ϕaeT,βyϕaa,ϕae=exp⁡j2πdλsin⁡ϕaesin⁡ϕaa,where ξ(t) is the noise vector, the elements of which are assumed to be stationary, zero-mean Gaussian random processes. And yk(t) is the waveform of a planewave arriving from the target DOA (ϕkaa,ϕkae). yk(t) may be specified as(3)ykt=εk∑j=1Mtajθkdujt,where uj(t) is a unit-norm waveform transmitted by the jth element and εk is the scatter (target) reflection coefficient. The effect of the multiple waveform transmission is given in aj(θkd) with j=1,…,Mt as the complex coefficient introduced by transmission of the jth waveform uj from the jth array element to kth target in the direction of departure (DOD) θkd.As shown in [22–25], it is often sufficient to consider that the mutual coupling coefficients between adjacent elements are almost the same and the magnitude of the mutual coupling coefficients decreases with increasing the antenna element spacing. In this paper, we assume the one antenna element is only affected by the coupling of the antenna elements, the spacing of which is within P×d, where P is a threshold. Then, we have(4)C=C1C2⋯CP0⋯0C2C1C2⋯CP⋯0⋮C2C1C2⋯CP0CP⋮C2⋱⋱⋯CP0CP⋮⋱C1C2⋮⋮0CP⋮C2C1C20⋯0CP⋯C2C1MrNr×MrNr,where C1,C2,…,CP are submatrices of C and are written as(5)C1=toeplitz1,c2d,c5d,…,cP×d,0,…,01×Nr,C2=toeplitzc2d,c5d,…,cP×d,0,…,01×Nr,⋮CP=toeplitzcP×d,0,0,…,01×Nr,where toeplitz{·} is the Toeplitz operator. And cγd is the mutual coupling coefficients between two antenna elements, the spacing of which is γd, where γ can be 2,5,…,P.In HF OTH radar, the application of DOA estimation or the adaptive beamforming is often in the domain of the two-dimensional range-Doppler spectra, where the “spike” peak of the target spectra is approximately occupied by several spectrum bins [26], one of which has maximum SNR within the target spectrum region. Therefore, we select the spectrum bin with the maximum SNR as the single snapshot of each of the receive antenna elements. And the antenna array vector x(t) is rewritten as x(0). Thus, the covariance matrix of the receiving signals can be written as [27](6)Rx=x0xH0=∑i=1KλieieiH+∑i=K+1Nr×MrλieieiH=EsΛsEsH+EnΛnEnH,where λ1>λ2>⋯>λK≥λK+1≥⋯≥λNr×Mr are the eigenvalues of Rx, e1,e2,…,eK,…,eNr×Mr are their associated eigenvectors, and Λs=diag⁡{λ1,…,λK} and Λn=diag⁡{λK+1,…,λNr×Mr}. The columns of Es and En are the eigenvectors associated with K largest eigenvalues and (Nr×Mr)-K smallest eigenvalues, respectively. For simplicity, we assume that the number of sources is known in this paper. With the known mutual coupling coefficients and the uncorrelated signals which have the high SNR, it is shown that the subspace spanned by Ca(ϕ1aa,ϕ1ae), Ca(ϕ2aa,ϕ2ae),…,Ca(ϕKaa,ϕKae) is approximately identical with the signal subspace spanned by Es and is approximately orthogonal to the noise subspace spanned by En. According to MUSIC algorithm, we have(7)EnHCaϕiaa,ϕiae≈0,i=1,…,K,PMUϕaa,ϕae=1EnHCaϕaa,ϕae2,where ϕaa and ϕae are the searching angles and EnHCa(ϕaa,ϕae) is zero whenever (ϕaa,ϕae)=(ϕiaa,ϕiae) and the peaks of PMU(ϕaa,ϕae) will correspond to the true DOAs.However, if the mutual coupling coefficients are unknown, the performance of DOA estimation would be degraded and the peaks ofPMU(ϕaa,ϕae) would point to the error DOAs. Moreover, the different signals can be represented as one of them multiplied by a complex factor in the case of single snapshot, which can be taken as coherent sources. Therefore, the mutual coupling correction technique and the spatial smoothing technique [28] are used to improve the performance of DOA estimation in the following part.Thex-axis and y-axis antenna elements of the receive antenna array are divided as the subarrays separately. Denote the data vector of the lth (l=1,2,…,Lx) subarray of the ith row along x-axis as xlxi(0), where Lx is the number of the subarrays, and the data vector of the nth (n=1,2,…,Ly) subarray of the jth column along y-axis as xnyj(0), where Ly is the number of the subarrays. We have(8)xlxi0=CAlxis0+ξlxi0,xnyj0=CAnyjs0+ξnyj0.And the covariance matrices are written as(9)Rfx=∑i=0Nr-1∑l=0Lx-1xlxi0xlxi0H,Rbx=∑i=0Nr-1∑l=0Lx-1Jxlxi0xlxi0HJ,Rfy=∑j=0Mr-1∑n=0Ly-1xnyj0xnyj0H,Rby=∑j=0Mr-1∑n=0Ly-1Jxnyj0xnyj0HJ,Rbfx=12Rfx+Rbx,Rbfy=12Rfy+Rby,where Rfx, Rbx, and Rbfx are the forward spatial smoothing covariance matrix, the backward spatial smoothing covariance matrix, and the forward-backward spatial smoothing covariance matrix of x-axis, respectively. Similarly, Rfy, Rby, and Rbfy are the forward spatial smoothing covariance matrix, the backward spatial smoothing covariance matrix, and the forward-backward spatial smoothing covariance matrix of y-axis. And J is the “exchange” matrix, which is denoted as(10)J=0…1⋮⋰⋮1⋯0.We denote the equivalent equations of the data vector after applying the spatial smoothing technique as(11)xbfx0=CbfxAbfxsbfx0+ξbfx0,xbfy0=CbfyAbfysbfy0+ξbfy0,where xbfx(0), Cbfx, Abfx, sbfx(0), and ξbfx(0) are the equivalent data vector, the equivalent mutual coupling coefficients matrix, the equivalent steering matrix, the equivalent signal vector, and the equivalent noise vector of x-axis, respectively. Similarly, xbfy(0), Cbfy, Abfy, sbfy(0), and ξbfy(0) are the equivalent data vector, the equivalent mutual coupling coefficients matrix, the equivalent steering matrix, the equivalent signal vector, and the equivalent noise vector of y-axis.In HF radar, the mutual coupling coefficients are unknown and the SNR is often not high. Therefore, the subspace spanned byCbfxabfx(ϕ1aa),Cbfxabfx(ϕ2aa),…,Cbfxabfx(ϕKaa) is not orthogonal to the noise subspace spanned by Enxbf, and the subspace spanned by Cbfyabfy(ϕ1ae),Cbfyabfy(ϕ2ae),…,Cbfyabfy(ϕKae) is not orthogonal to the noise subspace spanned by Enybf. Then, we have(12)EnxbfCbfxabfxϕ1aa,abfxϕ2aa,…,abfxϕKaa=0+ςx,EnybfCbfyabfyϕ1ae,abfyϕ2ae,…,abfyϕKae=0+ςy,where ςx and ςy are the residual items. Let # denote the pseudoinversion, and then we can get the least square solution as(13)Cbfx=Enxbf#ςxabfxϕ1aa,abfxϕ2aa,…,abfxϕKaa#,Cbfy=Enybf#ςyabfyϕ1ae,abfyϕ2ae,…,abfyϕKae#.By using the original DOA estimates ϕ-1aa,ϕ-2aa,…,ϕ-Kaa and ϕ-1ae,ϕ-2ae,…,ϕ-Kae obtained through the MUSIC algorithm, we can estimate the mutual coupling coefficients matrix by (13).


## 3. An Online Particle Mean-Shift Approach

We have assumed that the additive noise is stationary, zero-mean Gaussian random processes. From (13), we know that the mutual coupling coefficients matrix can be obtained from the DOAs we estimated. Thus, each of the elements of the estimated mutual coupling coefficients matrix is a random process. Therefore, we can now formulate a tracking problem in a Bayesian context for the online estimation for the mutual coupling coefficients matrix. In this section, we complete the model by introducing the idea of the particle filter and the mean-shift approach.From [29, 30], we know the residual items ςx and ςy in (13) can be represented as the linear combination of the additive noise. Thus, they are also stationary, zero-mean Gaussian random processes. And the bias of the estimates ϕ-1aa,ϕ-2aa,…,ϕ-Kaa and ϕ-1ae,ϕ-2ae,…,ϕ-Kae is Gaussian distribution [29, 30]. Denote Mrs and Nrs as the number of the elements of the subarrays on the x-axis and y-axis separately. We wish to calculate the posterior probability distributions p(Cbfx∣xbfx(0),ςx,ϕ-1aa,ϕ-2aa,…,ϕ-Kaa) and p(Cbfy∣xbfy(0),ςy,ϕ-1ae,ϕ-2ae,…,ϕ-Kae) which need the Mrs+K+1 and Nrs+K+1 dimension particles according to the idea of the particle filter algorithm [31] separately. To reduce the computational cost, (11) is rewritten as(14)xbfx0=C-bfxAbfxsbfx0+ΔC-bfxAbfxsbfx0+ξbfx0,xbfy0=C-bfyAbfysbfy0+ΔC-bfyAbfysbfy0+ξbfy0,where C-bfx and C-bfy are Toeplitz matrices and ΔC-bfx and ΔC-bfy are the residual items which satisfy(15)Cbfx=C-bfx+ΔC-bfx,Cbfy=C-bfy+ΔC-bfy.For any ϕaa and ϕae, we have(16)C-bfxabfxϕaa=Tabfxϕaac-x,C-bfyabfyϕae=Tabfyϕaec-y,where (17)c-x=1c1x⋯cMrs-1xT,c-y=1c1y⋯cNrs-1yT.Substitute (15) into (16), and we have (18)Cbfxabfxϕaa=Tabfxϕaac-x+ΔC-bfxabfxϕaa,Cbfyabfyϕae=Tabfyϕaec-y+ΔC-bfyabfyϕae.For each of the two steering vectors, such as abfx(ϕaa), T[abfx(ϕaa)] is a Mrs×Mrs matrix, which can be expressed as follows [22]:(19)Tabfxϕaa=T1abfxϕaa+T2abfxϕaa,T1abfxϕaai,j=abfxϕaai+j-1,i+j≤Mrs+1,0,otherwise,T2abfxϕaai,j=abfxϕaai-j+1,i≥j≥2,0,otherwise.According to (19), T[abfy(ϕae)] is also easily calculated, due to the fact that(20)EnxbfHCbfxabfxϕkaa=EnxbfHTabfxϕ-kaac-x+EnxbfHTabfxΔϕkaac-x+EnxbfHΔC-bfxabfxϕ-kaa+Δϕkaa+EnxbfHξbfx0fkx≈0,EnybfHCbfyabfyϕkae=EnybfHTabfyϕ-kaec-y+EnybfHTabfyΔϕkaec-y+EnybfHΔC-bfyabfyϕ-kae+Δϕ-kae+EnybfHξbfy0fky≈0,where Δϕkaa, Δϕkae, ΔC-bfx, and ΔC-bfy are the residual items with k=1,2,…,K. And EnxbfHξbfx(0)fkx and EnybfHξbfy(0)fky are the noise items, where fkx and fky are the factors. Define the coefficient matrices of the equations as(21)Qx=EnxbfHTabfxϕ-1aa⋮EnxbfHTabfxϕ-Kaa,Qy=EnybfHTabfyϕ-1ae⋮EnybfHTabfyϕ-Kae.Then, we have(22)Qxc-x=ςTx,Qyc-y=ςTy,where ςTx and ςTy are the residual vectors, the kth element of which is denoted as(23)ςTxk=-EnxbfHTabfxΔϕkaac-x-EnxbfHΔC-bfxabfxϕ-kaa+Δϕkaa-EnxbfHξbfx0fkx,ςTyk=-EnybfHTabfyΔϕkaec-y-EnybfHΔC-bfyabfyϕ-kae+Δϕ-kae-EnybfHξbfy0fky.As the additive noise is a stationary, zero-mean Gaussian random process, ςTxk and ςTyk are also Gaussian distribution, where k=1,2,…,K. Denote(24)Qx=q1xq2x⋯qMrsx,Qy=q1yq2y⋯qNrsy.Note that the first elements of the vectors c-x and c-y are equal to 1, and (22) can be rewritten as(25)Qx:,2:Mrsc-x2:Mrs=-q1x+ςTx,Qy:,2:Nrsc-y2:Nrs=-q1y+ςTy,where(26)Qx:,2:Mrs=q2x⋯qMrsx,Qy:,2:Nrs=q2y⋯qNrsy,c-x2:Mrs=c1xc2x⋯cMrsxT,c-y2:Nrs=c1yc2y⋯cNrsyT.Thus, the least square solution is easily obtained:(27)c-x2:Mrs=-Qx:,2:Mrs#q1x+ςTx,c-x1=1,c-y2:Nrs=-Qy:,2:Nrs#q2y+ςTy,c-y1=1.In (27), we have reduced the computational cost from Mrs+K+1 and Nrs+K+1 dimensions to Mrs+1 and Nrs+1 dimensions. According to (27) and (20), we can use the idea of the particle filter to obtain the estimates of the mutual coupling coefficients matrix.To reduce the computational cost further, (23) is approximately equal to(28)ςTxk≈EnxbfHTabfxϕ-kaa+Δϕ~kaac-x-EnxbfHTabfxϕ-kaac-x,ςTyk≈EnybfHTabfyϕ-kae+Δϕ~kaec-y-EnybfHTabfyϕ-kaec-y,where Δϕ~kaa and Δϕ~kae are the residual items. And substituting (28) into (20) and (22), we have(29)EnxbfHCbfxabfxϕkaa≈EnxbfHTabfxϕ-kaa+Δϕ~kaac-x,EnybfHCbfyabfyϕkae≈EnybfHTabfyϕ-kae+Δϕ~kaec-y.Moreover, (29) can also be rewritten as(30)EnxbfHCbfxabfxϕkaa≈E^nxbfHTabfxϕ-kaa+Δϕ^kaac-x,EnybfHCbfyabfyϕkae≈E^nybfHTabfyϕ-kae+Δϕ^kaec-y,where E^nxbf and E^nybf are the noise space of the covariance matrices toeplitz{Rbfx} and toeplitz{Rbfy} separately, which are the Toeplitz matrix form of Rbfx and Rbfy by the Toeplitz approximation transformation [32]. According to (24)–(27), we can also obtain the estimate of the mutual coupling coefficients matrix by (29) or (30).For convenience, we denote the probability density functions of theith estimation by the idea of the particle filter as p(Cbfix∣ςTix) and p(Cbfiy∣ςTiy). That is to say, we define(31)pCbfx∣xbfx0,ςx,ϕ-1aa,ϕ-2aa,…,ϕ-Kaa⟶pCbfix∣ςTix,pCbfy∣xbfy0,ςy,ϕ-1ae,ϕ-2ae,…,ϕ-Kae⟶pCbfiy∣ςTiy.And p(Cbfix∣ςTix) and p(Cbfiy∣ςTiy) are the a posteriori probability density functions in Mrs×Mrs and Nrs×Nrs dimensionality space. To find a robust result of Cbfx and Cbfy, a mean-shift algorithm based on a kernel function is introduced in the following, which is a nonparametric density gradient estimation in nature. Conditions on the kernel function are derived to guarantee asymptotic unbiasedness, consistency, and uniform consistency of the estimates. This approach is a simple iterative procedure that shifts each data point to the average of data points around its neighborhood. And the mean-shift approach is applied to search the most rapid increase of the probability density function for the matrix tracking in the case of the estimated mutual coupling coefficients matrices generated by the particles.Denote a new kernel function estimation ofp(Cbfix∣ςTix) and p(Cbfiy∣ςTiy) [32, 33], and we have(32)pFCbfx∣ςTx=∑i=1LFCbfix-Cbfx/hxpCbfix∣ςTixhxMrs×Mrs∑i=1LpCbfix∣ςTix,pFCbfy∣ςTy=∑i=1LFCbfiy-Cbfy/hypCbfiy∣ςTiyhyNrs×Nrs∑i=1LpCbfiy∣ςTiy,where L is the number of the particles and hx and hy are the radius of the high-dimensional ball areas in the Mrs×Mrs- and Nrs×Nrs-dimensional Euclidean space separately. And the Epanechnikov kernel can be expressed as follows [33, 34]:(33)FCbfx=1-Cbfx2,ifCbfx≤1,0,ifCbfx>1,FCbfy=1-Cbfy2,ifCbfy≤1,0,ifCbfy>1.Let (34)GCbfx=-F′Cbfx,GCbfy=-F′Cbfy.The gradient estimation of p(Cbfx∣ςTx) and p(Cbfy∣ςTy) is(35)∇^pCbfx∣ςTx=∇p^FCbfx∣ςTx=2∑i=1LCbfx-CbfixF′Cbfix-Cbfx/hx2ωCbfixhxMrs×Mrs+2∑i=1LωCbfix,∇^pCbfy∣ςTy=∇p^FCbfy∣ςTy=2∑i=1LCbfy-CbfiyF′Cbfiy-Cbfy/hy2ωCbfiyhyNrs×Nrs+2∑i=1LωCbfiy,where ω(·) is the weight function. Equation (35) can be rewritten as(36)∇^pCbfx∣ςTx=2∑i=1LCbfix-CbfxGCbfix-Cbfx/hxωCbfixhxMrs×Mrs+2∑i=1LωCbfix=2hx2∑i=1LGCbfix-Cbfx/hxωCbfixhxMrs×Mrs+2∑i=1LωCbfix∑i=1LCbfix-Cbfx/hxGCbfix-Cbfx/hxωCbfix∑i=1LGCbfix-Cbfx/hxωCbfix,∇^pCbfy∣ςTy=2∑i=1LCbfiy-CbfyGCbfiy-Cbfy/hyωCbfiyhyNrs×Nrs+2∑i=1LωCbfiy=2hy2∑i=1LGCbfiy-Cbfy/hyωCbfiyhyNrs×Nrs+2∑i=1LωCbfiy∑i=1LCbfiy-Cbfy/hyGCbfiy-Cbfy/hyωCbfiy∑i=1LGCbfiy-Cbfy/hyωCbfiy.Denote(37)p^GCbfx∣ςTx=∑i=1LGCbfix-Cbfx/hxωCbfixhxMrs×Mrs+2∑i=1LωCbfix,p^GCbfy∣ςTy=∑i=1LGCbfiy-Cbfy/hyωCbfiyhyNrs×Nrs+2∑i=1LωCbfiy.Then, we use Mh(Cbfx) and Mh(Cbfy) to express the formula followed p^G(Cbfx∣ςTx) and p^G(Cbfy∣ςTy), and we have(38)MhCbfx=∑i=1LCbfix-Cbfx/hxGCbfix-Cbfx/hxωCbfix∑i=1LGCbfix-Cbfx/hxωCbfix,MhCbfy=∑i=1LCbfiy-Cbfy/hyGCbfiy-Cbfy/hyωCbfiy∑i=1LGCbfiy-Cbfy/hyωCbfiy,where p^G(Cbfx∣ςTx) and p^G(Cbfy∣ςTy) are the estimations of kernel function G(C) for p(Cbfx∣ςTx) and p(Cbfy∣ςTy) and Mh(Cbfx) and Mh(Cbfy) are the mean-shift vectors. Then, (35) could be rewritten as(39)∇^pCbfx∣ςTx=∇p^FCbfx∣ςTx=2hx2p^GCbfx∣ςTxMhCbfx,∇^pCbfy∣ςTy=∇p^FCbfy∣ςTy=2hy2p^GCbfy∣ςTyMhCbfy,(40)MhCbfx=2hx2∇p^FCbfx∣ςTxp^FCbfx∣ςTx,MhCbfy=2hy2∇p^FCbfy∣ςTyp^FCbfy∣ςTy.We can conclude from (40) that the vectors Mh(Cbfx) and Mh(Cbfy) are proportional to the gradient of p^F(Cbfx∣ςTx) and p^F(Cbfy∣ςTy). Hence, Mh(Cbfx) and Mh(Cbfy) point in the direction of the most rapid increase of the probability density function.The routine of the online particle mean-shift approach is shown in Algorithm1. We refer to this algorithm as OPMA for convenience.Algorithm 1 (OPMA).
InitializeCbfx and Cbfy as the identity matrices; this initial DOA’s estimates are calculated by the MUSIC algorithm with the covariance matrices obtained by (9). And the sums of the MUSIC spectrum amplitude of the peaks of the targets are denoted as Px(0) and Py(0).
Fori=0,1,2,…,L, do (1)
LetΔϕ~kaa and Δϕ~kae in (29) and Δϕ^kaa and Δϕ^kae in (30) be the random particles generated by rand(1)·Δθ~ (k=1,2,…,K), where rand(1) is a stochastic number with the uniform distribution within the range of -1 to 1 and Δθ~ is a quarter of the corresponding beamwidth.(2)
Substitute the random particles into (29) and (30), and the coupling coefficient matrix estimates Cbfx and Cbfy are calculated by (27). Then, the MUSIC spectra are obtained by the noise spaces Enxbf, Enybf, E^nxbf, and E^nybf and the coupling coefficient matrix estimates. And the sums of the MUSIC spectrum amplitude of the peaks of the targets are denoted as Px(i), Py(i), P^x(i), and P^y(i).(3)
IfPx(0)>Pthreshold, let SPx(i)=Px(i). Otherwise, let SPx(i)=P^x(i), where SPx(i) is a temporary variable. If SPx(i)>Px(0), we reserve the results SPx(i) for the calculation of the mean-shift algorithm. Otherwise, the results are abandoned. Similarly, if Py(0)>Pthreshold, let SPy(i)=Py(i). Otherwise, let SPy(i)=P^y(i), where SPy(i) is also a temporary variable. If SPy(i)>Py(0), we reserve the results SPy(i) for the calculation of the mean-shift algorithm. Otherwise, the results are abandoned
end for(4)
The mean-shift results are calculated with the results selected in step (3) by (38). And we select the weights ω(Cbfix) and ω(Cbfiy) in (38) as SPx(i)/∑i=0LSPx(i) and SPy(i)/∑i=0LSPy(i).(5)
Select the mean-shift results as the estimates of the equivalent mutual coupling coefficients matrix ofx-axis and y-axis. Then, substitute into the MUSIC algorithm separately, and the estimates of the azimuth and the elevation of the targets can be obtained.Remark 2.
P t h r e s h o l d is an empirical value, which is the spatial spectrum amplitude of the MUSIC algorithm and can be about 10 dB. And the number of the particles L can be about 20, which has the moderate computational complexity.


## 4. Simulations

In this section, we evaluate the performance of the proposed approach mentioned above. Consider the receiver antenna array of a MIMO HF sky-wave radar system, which is a 2D uniform rectangle array with 50 identical and omnidirectional antennas inx-axis and y-axis separately, and the ratio of the distance between the neighboring elements d and the wavelength λ is 0.5. Assume p=13, and c1d=0.7+0.8j, c2d=0.5+0.6j, c3d=0.1+0.2, c4d=0.01+0.02j, c2d=0.6+0.7j, c5d=0.4+0.5j, c8d=0.2+0.3j, c10d=0.07+0.08j, and c13d=0.03+0.04j. And the receiving targets have the same range which is 1000 km, and they have the same DOD and the same radial velocity. The snapshot is single one and the noise is the additive white Gaussian noise. Four approaches are considered in the following simulations, which are the spatial smoothing technique without the mutual coupling correction (which is referred to as SST), the spatial smoothing technique with the mutual coupling correction by (25) with ςTx=0 and ςTy=0 (which is the Conventional Mutual Coupling Coefficients Estimation Method and is referred to as CMCCE), the spatial smoothing technique with no mutual coupling, and the OPMA proposed in this paper. We consider the success probability, the estimation bias, and the estimation variance of the four approaches in the case of SNR = −10 dB~30 dB. 100 simulations are carried out for each different SNR. The estimation bias in each Monte Carlo simulation is measured by the root mean squared error (RMSE).In the first simulation, we consider three targets with the different azimuths which are −4°, −1°, and 15° and the same elevation which is 90°. The success probability of the first simulation is shown in Figure2. It is shown that OPMA has higher success probability than SST and CMCCE and lower success probability than the no coupling result in the low SNR region. The RMSE of the first simulation is shown in Figure 3. It is shown that OPMA has lower RMSE than SST and CMCCE and higher RMSE than the no coupling result. The variance of the first simulation is also calculated. The result is shown in Figure 4. It is shown that the proposed approach has smaller variance than SST and CMCCE.Figure 2
Success probability of the first simulation.Figure 3
RMSE of the first simulation.Figure 4
Variance of the first simulation.In the second simulation, we consider three targets with the different elevations which are 0°, 3°, and 18° and the same azimuth which is 0°. The success probability of the second simulation is shown in Figure5. It is shown that OPMA has higher success probability than SST and CMCCE and lower success probability than the no coupling result in the low SNR region. The RMSE of the second simulation is shown in Figure 6. It is shown that OPMA has lower RMSE than SST and CMCCE and higher RMSE than the no coupling result. The variance of the second simulation is also calculated. The result is shown in Figure 7. It is shown that the proposed approach has smaller variance than SST and CMCCE.Figure 5
Success probability of the second simulation.Figure 6
RMSE of the second simulation.Figure 7
Variance of the second simulation.


## 5. Conclusion

We introduce the OPMA for the mutual coupling coefficients estimation of the 2D uniform rectangle array with the single snapshot in MIMO HF sky-wave radar. The equation of the mutual coupling coefficients matrix and the spatial smoothing technique of the uniform rectangle array are introduced. Besides, the equivalent equations of the mutual coupling coefficients matrix after the preprocessing of the spatial smoothing technique are derived. Based on the idea of the particle filter and the mean-shift algorithm, the OPMA is proposed, which has a moderate computational complexity. We compare the performance of the proposed approach with the no mutual coupling case and the approaches of SST and SMCCE. The proposed approach is more robust and more accurate than SST and SMCCE. The simulation results demonstrate the validity of the proposed method in terms of the success probability, the statistics of bias, and the variance.



---

*Source: 7515028-2016-07-20.xml*