A new algorithm to estimate the direction of arrival (DOA) and polarization parameters of signals impinging on an array with electromagnetic (EM) vector-sensors is presented by exploiting the canonical polyadic decomposition (CPD) of tensors. In addition to spatial and temporal diversities, further information from the polarization domain is considered and used in this paper. Estimation errors of these parameters are evaluated by the Cramér-Rao lower bound (CRB) benchmark, in the presence of additive white Gaussian noise (AWGN). The superiority of the proposed algorithm is shown by comparing with the derivative algorithms of MUSIC and ESPRIT. In the proposed algorithm, the parameters can be estimated by virtue of the diversities of the spatial and polarization belonging to the factor matrices, rather than the conventional subspace which is the foundation of MUSIC and ESPRIT. Additionally, the classical CPD algorithm based on Alternating Least Squares (ALS) is introduced to verify the efficacy of the proposed CPD algorithm. Results demonstrate that when the number of snapshots is greater than 50, the proposed algorithm requires a smaller number of snapshots to achieve a high level of performance, compared against the subspace-based algorithms and the ALS-based algorithm. Furthermore, in the matter of the array with a small number of sensors, the discovered advantage concerning the Root Mean Square Error (RMSE) in estimating the DOA and the polarization state of the signal is noteworthy.

1. Introduction

Signal parameters estimation is of great significance in many applications such as satellite navigation, wireless communication, radar, and sonar. An array of multiple sensors placed in different spatial locations is used to estimate the signal parameters, which mainly include spatial and polarization parameters of the signal. The spatial parameter refers to the direction of arrival (DOA) of the signal, while the polarization auxiliary angle and the polarization phase difference of the signal can be regarded as the polarization parameters. Among signal parameters estimation algorithms based on matrix operations, some depend on exhaustive search, e.g., multiple signal classification (MUSIC) [1], whereas others do not, e.g., estimation of the signal parameters via rotational invariance techniques (ESPRIT) [2] and root-MUSIC [3]. One of the common features of the matrix-based algorithms is that the matrix can only reflect the two-dimensional diversity of the desired signal. Some methods exploiting the signal-subspace embodied by MUSIC, ESPRIT and their derivative algorithms have been introduced in [412], in which the tensor-based data model has been combined in [1012]. Based on the Higher-Order Singular Value Decomposition (HOSVD) of the measurement tensor, two algorithms of parameters estimation combined with ESPRIT have been developed in [10]. Similarly, in [11, 12], several tensor MUSIC methods based on HOSVD have been derived. The existing research results show that the performance of the subspace-based tensor method is improved compared to the subspace-based matrix method, and the computational complexity is higher in the tensor method than in the matrix method but of the same order. In [13], a greedy algorithm called randomized multiple candidate iterative hard thresholding for DOA estimation has been proposed. A set of potential candidates using the iterative hard thresholding algorithm are first obtained, and then the best candidate can be selected based on the a priori knowledge of the distribution of the signal and noise matrices.

The canonical polyadic decomposition (CPD) of the third-order tensor is a minimal decomposition into a sum of rank-1 tensors [14]. The uniqueness of CPD, which is often called essential uniqueness in engineering papers, has been extensively studied in the field of algebraic geometry. The property of uniqueness makes the CPD a basic tool for signal separation, and it is widely used in telecommunication, array processing, machine learning, etc. [1520]. In [15], a novel multiuser receiver for joint symbol and channel estimation by capitalizing on a tensor modeling of the end-to-end system has been proposed. The multiple invariance sensor array processing has been linked to the uniqueness of CPD in [20]. It shows the uniqueness of single and multiple invariance ESPRIT, which stems from a deterministic decomposition of the third-order tensor signal model. Moreover, it provides the DOA of the source with less restrictions of signal stationarity than the aforementioned statistical approaches. In [21], a method of DOA estimation for seismic plane waves has been considered from a deterministic perspective using CPD. In addition to temporal and spatial information, the different propagation speed of waves is taken into account by using the multidimensional feature of the tensor. In [22], a CPD-based approach for distinguishing the signals with the same DOA and copolarized state has been proposed. However, the computational complexity of the algorithm is very high due to the complexity of the problem model.

The main purpose of this paper is to explore deterministic decomposition of the third-order tensor signal model by exploiting the uniqueness of CPD and then estimate the signal parameters. The polarization diversity of the signal is considered in this paper, in addition to the temporal and spatial sampling employed by the most current parameters estimation methods. Note that the proposed algorithm focus on not only the DOA estimation but also the estimation of polarization parameters delivered by the signal. The strength of our approach is its excellent performance compared against the traditional subspace-based algorithms such as MUSIC, HOSVD-MUSIC, and HOSVD-ESPRIT given shorter snapshots and the array with a small number of sensors.

The rest of this paper is organized as follows. Section 2 reviews the CPD prerequisites. The tensor signal model is derived in Section 3, and the introduction to traditional subspace-based algorithms is also provided. The CPD-based algorithm is proposed in Section 4, followed by the numerical simulations in Section 5. Finally, the last part presents our concluding remarks. The algebra notations involved in the paper are shown in Table 1.

2. Prerequisites for CPD

Let the -th entry of the third-order tensor be symbolized by . Herein, we define ; i.e., there exist three nonzero vectors , , and such that , which means that stands for all possible values of indices.

A third-order tensor with general-rank can be expressed as a sum of rank-one items:where , , and . The decomposition shown in (1) is called the polyadic decomposition (PD). Particularly, if the number of rank-one terms is minimal, then is the rank of (i.e., ), and (1) is defined as the CPD of .

Additionally, (1) can also be written as , where the matrices , , and are defined as the first, second, and third factor matrix of , respectively. Obviously, the entries of can be written as . In addition, here is a potential scenario: the CPD of a general-rank tensor is not unique. Nevertheless, we still call that the CPD of is unique when it is only subject to the trivial indeterminacies, which include that the rank-one terms shown in (1) can be arbitrarily permuted and the vectors belonging to any rank-one term can be arbitrarily scaled. Similarly, the factor matrices for any two CPDs and coincide up to column permutation and scaling; for instance, can be obtained by scaling and permuting the columns of (see [9] for more details). Herein, for reshaping the tensor into a matrix , we introduce the following definition.

Definition 1. The matrix unfolding that yields and the permutation , for the tensor is defined as . The -th entry of is given bywhere for , and , . Assuming that and the permutation , the corresponding matrix unfolding of (1) can be expressed aswhere , symbolizes the th horizontal slice of , denotes the th row of , and is the diagonalization operator. What is more, in the case of , we can obtain the other matrix unfoldings:where and . This paper mainly uses the algebraic algorithms described in [14, 23] to explore the application of CPD in signal parameters estimation. The key of the algebraic algorithms presented in [14, 21] is the following theorem, which has been derived in [24].

Theorem 2. Define a third-order tensor , and assume that the factor matrices and have full column rank and any two columns selected from the factor matrix are linearly independent:where represents the largest number such that every subset with columns of the factor matrix is linearly independent. Then ; the CPD of can be computed by algebraic algorithms uniquely.

One algebraic algorithm elaborated in [23] is adopted in this paper. Herein, we consider the particular case of that conforms to Theorem 2. By (3), the matrix unfolding of can be expressed aswhich obtains the equations that

Since , the diagonal elements of the diagonal matrix are distinct. Therefore, the factor matrices and can be uniquely recovered by the eigenvalue decomposition (EVD) of (7), and then the matrix can be easily identified from (6). However, the factor matrices are not necessarily square matrices in reality. Hence, we need to consider the generalized inverse operation of the matrices to obtain the CPD of the tensor.

Herein, in order to avoid ambiguity in signal modeling, some of the main assumptions are stated as follows.

A1. Completely polarized wave: the horizontal and vertical electric components located in the wavefront of the transverse electromagnetic (TEM) wave are completely correlated.

A2. Far-field assumption: we assume that the distance between the source and the receiving array is much larger than the array aperture. Thus, the source of the signal can be considered point-like. Moreover, the wavefront of the TEM wave can be approximated as a plane at the sensor level.

A3. Narrow-band assumption: the reciprocal of the signal bandwidth is much greater than the maximum propagation time of the signal passing through the array aperture. Typically, the signal received by an array is a varying complex envelope, which is modulated at a carrier with high frequency (HF). In order to facilitate the digital signal processing, the intermediate frequency (IF) or baseband signal is usually wanted. Hence, the analog down-conversion for the received HF signal is a more common practice.

3. Signal Modeling

This section mainly focuses on the three-dimensional feature, i.e., the vector-sensor enjoys in spatial, polarization, and time diversities and then establishes the third-order tensor-based signal model. In addition, conventional matrix-based solutions are introduced briefly. Firstly, the diversities involved in this paper are shown as follows.

Time Diversity. The baseband analog signal obtained by down-conversion is a time-dependent function. Moreover, the discrete-time signal can be obtained by sampling the baseband analog signal in time domain. It is assumed that signals with the complex amplitudes impinge on the vector-sensor array. Thus, the time diversity can be formulized aswhere , represents the snapshot at the instant .

Spatial Diversity. The spatial sampling of the array reflects the spatial diversity. Consider a vector-sensor array composed of EM vector-sensors located at for . Combined with the time diversity, the signal received by the array can be expressed aswhere the symbol denotes the wave length of the signal impinged on the vector-sensor array. Moreover, denotes the incident direction of the th signal source associated with the azimuth angle and the elevation angle .

Polarization Diversity. The EM vector-sensor is usually composed of multiple magnetic and electric short dipoles in a common spatial point. Hence, the EM vector-sensor is capable of measuring multiple electromagnetic components with the number of not less than one. For example, a three-dipole sensor can measure the three components of the electric field in three perpendicular directions. Furthermore, the complete EM vector-sensor has the ability to obtain all electromagnetic components (i.e., three electric components and three magnetic components). In particular, if the vector-sensor array is composed of a monopole, the orientations of the monopoles cannot be the same; otherwise, the array will degrade into a scalar-array. We formulate the polarization diversity combined with the time diversity,where , is the number of component measured by an EM vector-sensor, for three-dipoles, and for complete EM vector-sensor. Moreover, represents the th component of the polarization steering vector , and the symbol denotes the spatial-polarization parameter of the th signal source with respect to the polarization auxiliary angle and the polarization phase difference .

3.1. Tensor-Based Signal Model

Consider the diversities mentioned above, and assume that all the signals are cofrequency and incoherent. The spatial steering vector of the signal with DOA is expressed as

The polarization steering vector of the signal with spatial-polarization parameter is given bywhere and represents the polarization selection matrix of the EM vector-sensor. Particularly, represents the complete EM vector-sensor, and symbolizes the three-dipole sensor, where denotes the identity matrix. Additionally, the first three rows of the matrix represent the projections of the electric field component of the incident signal with DOA in the three axial directions of the coordinate system, and the last three rows of the matrix represent the projections of the magnetic field component in the three axial directions. The spatial-polarization steering vector can be further defined as

Hereby, let denote the spatial white noise at the instant , which is additive and Gaussian complex circular. The vector-output of the array at the instant can be expressed as

Considering that the multidimensional nature of the tensor just caters to the multidomain diversities of the signal, the following spatial-polarization steering matrix can be obtained by performing the outer product operations of (12) and (11):where . Let denote the consecutive snapshots of the signal; the tensor-output of snapshots can be modeled aswhere . The tensor noise is expressed as , which is yielded by the tensorization with respect to the noise matrix .

3.2. Traditional Subspace-Based Algorithms

Subspace decomposition methods with respect to EVD and SVD are the mainstream of MUSIC, ESPRIT, and their derivative algorithms. Herein, we introduce the matrix-based MUSIC algorithm and the tensor-based MUSIC algorithm based on HOSVD. The covariance matrix for snapshots corresponding to (14) can be given bywhere . Then the eigenvector corresponding to the small eigenvalues of can be written in a matrix form , where . Moreover, according to (11) and (12), we need to construct a search matrix:

The solution of the matrix-based MUSIC algorithm can be defined aswhere . The estimation of the polarization parameter is detailed in Appendix A. For the tensor-based MUSIC algorithm based on HOSVD [11], the following signal model needs to be constructed:where . The covariance tensor for snapshots is expressed aswhere . Then, is decomposed using the HOSVD procedure aswhere represents -mode product, is the core tensor, and can be obtained by SVD of -mode matrix unfolding of tensor . Finally the parameters of the sources are obtained by maximizing the following criterion:where is a truncated form of (see [11] for specific truncation methods). In addition, the HOSVD-ESPRIT algorithm depending on translation invariant property of the array is elaborated in Appendix B.

4. CPD-Based Algorithm

We hereby execute the CPD of the tensor-based signal model to obtain the factor matrices, which are then used to extract the parameters of the signals. Let , where , , and . Assume that at least two factor matrices of have full column rank and the -rank of the rest factor matrix is not less than 2, which are consistent with the Theorem 2. Moreover, there are conditions in which , , and are not less than , respectively, and this cannot be ignored here to ensure the uniqueness of the factor matrices.

4.1. Computation of CPD

As is well known, CPD of a third-order tensor can be achieved by the classical Alternating Least Squares (ALS) algorithm, which is detailed in Appendix C. Differently, we hereby propose an EVD-based algebraic method to achieve CPD of tensor. The tensor-based signal model (16) can be unfolded aswhere denotes the th row of the , . According to (24), the following relationship can be easily derived:where , , and denotes the Moore-Penrose inverse. From (24), we have , and . Hence, for any , . The features of the CPD are demonstrated as follows:(1)The columns of the factor matrix are the eigenvectors of . The corresponding eigenvalues are the first larger eigenvalues, which are arranged in a descending order.(2)The eigenvectors corresponding to the -large eigenvalues of constitute the factor matrix .(3)The rows of factor matrix can be recovered thanks to the relationship .

Obviously, when , there are various combinations that can be used for EVD. In order to make full use of the statistical properties of (24), we select possible combinations for EVD to obtain , , where . We follow the convention that the -tuples are ordered lexicographically: the th tuple is preceding the th tuple if and only if , . According to the trivial indeterminacies aforementioned for CPD uniqueness, we know that and , , coincide up to column permutation and scaling. In order to achieve the addition of multiple sets of results, we need to match the multiple sets of results. That is to say, the eigenvectors in each factor matrix obtained by the EVD need to be arranged in a descending order according to the corresponding eigenvalues, and all the eigenvectors need to be normalized. After the matching process, and can be obtained, and then the following statistical operations are performed:

Next, the factor matrix can be recovered according to the results of (26).

4.2. Parameters Estimation from the Factor Matrices

From the aforementioned tensor-based signal model, we can see that the spatial-polarization parameters of signals can be reflected by the factor matrices and the factor matrix that only contains the spatial parameters of signals. Hence, once the factor matrices are obtained, we first need to extract the spatial parameters from .

4.2.1. Estimation of Spatial Parameters

Spatial DOA of the th signal impinging on the array can be calculated directly by solving the following equations:

Obviously, any two equations belonging to (27) can be used to derive the DOA of the th signal. Therefore, we can get sets of solutions from (27). Improve the accuracy of the estimated parameters by averaging multiple sets of solutions:where represents the solution of the th subequations . We follow the convention that the -tuples are ordered lexicographically: the th tuple is preceding the th tuple if and only if , . Furthermore, we can construct the optimization problem as follows to estimate the DOAs of the signals:where represents the th column of the factor matrix obtained by CPD, represents the amplitude vector of , and the operator denotes the 2-norm. The solution method to (29) and its performance can be used as a further work, and this paper discusses only the performance based on the solution of (28).

4.2.2. Estimation of Polarization Parameters

Once the estimation of DOA is obtained, the polarization parameters of the th signal can be extracted fromwhere represents the th column of the factor matrix obtained by CPD. Hereby, the estimated value of can be expressed as , and the estimated values of the polarization parameters can be obtained bywhere is the arctangent operator and the operator denotes phase acquisition.

Obviously, the estimation accuracy of the DOA only depends on the factor matrix obtained by CPD, and the estimation accuracy of the polarization parameters is affected by both the estimation accuracy of the DOA and the matrix obtained by CPD. It is well known that, in order to satisfy the validity of the parameters estimation, the algorithm must follow the premise of . According to [14], there exists Kruskal bound ensuring the CPD uniqueness . With the noise present, the accuracy of the CPD is severely impaired, when is very close to the bound. However, as increases, the dependence of the CPD accuracy on will rapidly weaken. In addition, as can be seen from (26) and (28), the proposed algorithm can obtain multiple sets of parameters estimation results in each CPD and perform statistical operations on these results. However, what is different is that the subspace-based methods requires sufficient statistical information, which is usually obtained by increasing the number of sensors or snapshots, to ensure the accuracy of the parameters estimation, even if the number of the parameters to be estimated is very small. We, therefore, conclude that the proposed algorithm can achieve relatively desirable performance by using a smaller number of snapshots as well as a smaller number of sensors compared with the subspace-based methods. The subsequent simulation results also justify the above theoretical analysis.

4.3. Computational Complexity

The parameters estimation algorithm using EVD-based CPD of tensors is summarized in Algorithm 1, which is easy to implement by row-by-row access. Note that the computational complexity scales linearly to the number of loops, since the procedure of the algorithm runs serially.

Require:  , , , , .
1: Reshape the matrcies by unfolding the tensor ;
2: For to do
3:Calculate , ;
4:Obtain by executing the EVD of ;
5:Obtain by executing the EVD of ;
6:Obtain and by matching processing.
7: End for
8: ;
9: Recover according to ;
10: For to do
11:For to do
12:Calculate by solving equations ;
13:End for
14: ;
17: End for
Ensure:  , .

In fact, the EVDs of (25) follow from the generalized EVDs (GEVDs) of the matrix pencils and . In order to evaluate the complexity of the proposed algorithm, we first need to know the number of multiplications required by GEVD. An effective GEVD solution is given by the method of iteration [25, pp.413-414], in which multiplications are required for two matrices with rank-, where the constant depends on the selected iterative strategy, e.g., for Hessenberg-Triangular Reduction, . The Moore-Penrose inverse of the matrix or matrix involves multiplications, where . Additional multiplications are also required to recover the factor matrix . Therefore, the computational complexity of the EVD-based CPD does not exceed , since the premise of has been established. Furthermore, in order to compare the computational complexity for the proposed algorithm, we analyze the required number of multiplications for subspace-based methods. For the aforementioned matrix-based method, the computational complexity of obtaining subspaces is approximately . Considering the operation process of HOSVD based on orthogonal iterations [25, pp.454-455], the computational complexity of subspace for tensor model is about . Obviously, the conclusion that the computational complexity of the CPD algorithm is lower than that of the tensor-based HOSVD methods can be obtained through the above analysis. Additionally, for the extraction process of the parameters, the computational complexity of the search strategies is generally much higher than that of the equation-solving methods, under the premise of ensuring the estimation accuracy. Therefore, the computational complexity of the proposed algorithm is lower than that of the search-based methods.

Additionally, according to the procedure of the ALS-based CPD described in Appendix C, the computational complexity of the factor matrix is not less than , since the premise of has been established. Similarly, the computational complexity of the factor matrices and is not less than . Therefore, the computational complexity of the ALS-based CPD is about , where represents the number of iterations. We use the Matlab 2016a to implement the after-mentioned simulations on a computer with Intel Core i5 CPU6200U 2.3GHz and 8GB memory running Windows 10. From the perspective of computing time, the time required for ALS-based CPD is about twice that of EVD-based CPD to obtain the subsequent simulation results.

5. Simulation Results

The superiorities of the proposed algorithm will be demonstrated by the results of numerical simulations in this section. We consider the following scenario with a uniform linear array (ULA) separated by half wavelength. The array consists of complete EM vector-sensors, and two signals with the polarization parameters and impinge on the array from and , respectively. Both signals satisfy the aforementioned assumptions A1 to A3. The modulation modes of the signals are BPSK and QBPSK, respectively, and the symbol transmission rate of both signals is . Additionally, modulation frequency and sampling frequency are 16kHz and 32kHz, respectively. The noise component is assumed to be zero-mean additive white Gaussian noise. Performance of the EVD-based CPD algorithm is compared with three other subspace-based algorithms, MUSIC, HOSVD-MUSIC, HOSVD-ESPRIT, and the ALS-based CPD algorithm. The performances of the algorithms are evaluated by the Root Mean Square Error (RMSE):where denotes one of the parameters , and is the estimation of in the th trial. Throughout all simulations, results are averaged by 500 Monte-Carlo trails and compared with the Cramér-Rao lower bound (CRB) benchmark, which is described in Appendix D.

Example 3. Assume that the number of complete EM vector-sensors is 4, and the signal-to-noise ratio (SNR) is 30dB. The proposed EVD-based CPD algorithm is compared with the subspace-based algorithms and the ALS-based CPD algorithm in terms of the DOA RMSE versus the different number of snapshots . Figure 1 illustrates a superior performance of the proposed EVD-based CPD algorithm with respect to the other algorithms when the number of snapshots is larger than 50. It is worth noting that, when the number of snapshots is greater than 100, the performance change of the ALS-based CPD algorithm with the increase of snapshots is negligible. In the case where the number of snapshots is extremely small, the main reason for the poor performance of the EVD-based CPD algorithm is that the number of the signals to be estimated is close to the uniqueness condition for CPD.

Example 4. Scenarios with various numbers of sensors () are considered in this example. Each of these scenarios is tested versus SNR, and the number of snapshots is fixed to . The following conclusions can be obtained by observing Figure 2.(1)For all algorithms in each case, RMSEs decrease gradually with the increase of SNR. Note that the EVD-based CPD algorithm performs better than the other algorithms from case to case when SNR is greater than dB.(2)In the case of low SNR, i.e., SNR is less than dB, the problem of the poor estimation accuracy for the proposed EVD-based CPD algorithm is effectively improved with the increase of the number of sensors.(3)The accuracy of the DOA estimation for MUSIC algorithm changes better obviously, as the number of sensors increases. This is due to the fact that the MUSIC algorithm determines the DOA by executing an exhaustive search [21].

Example 5. Once the estimation of DOA is obtained in Example 4, the polarization parameters can be estimated. Figure 3 shows the estimation accuracy of the parameter by the EVD-based CPD and the other algorithms in the scenarios in Example 4. The differences in performance between algorithms are similar to the tendencies emerged in Figure 2. Similarly, the performance of the ALS-based CPD is significantly improved as the number of sensors increases. However, the RMSEs of the algorithms for are far away from the CRB boundary, especially for and . This loss can be explained by the 2-step method of the estimation process, which consists of estimating DOA parameters firstly, and then corresponding polarization parameters in a second step. In other words, the accumulation of errors in the 2-step calculation leads to this loss.

6. Conclusion

In this paper, we proposed an EVD-based CPD approach for tensors to tackle the problem of signal parameters estimation. The procedure of the EVD-based CPD was described, followed by expatiations and simulations of the EVD-based CPD algorithm applied to the signal parameters estimation. The traditional subspace-based algorithms, MUSIC, HOSVD-MUSIC HOSVD-ESPRIT, and the classical ALS-based CPD algorithm, were introduced to verify the efficacy of the proposed EVD-based CPD algorithm.

The strength of the proposed algorithm lies in the estimation of the factor matrices from the tensor model of the signal, which is then used to extract the signal parameters. The signal parameters can be estimated by virtue of the diversities of the spatial and polarization belonging to the factor matrixes. Furthermore, when the number of sensors is small in the array, the proposed algorithm achieves a lower RMSE.


A. Polarization Parameters Estimation of MUSIC

Through (18) and (19), we can construct a new matrix:

The polarization parameters estimation of the signals can be expressed aswhere , , and the operator symbolizes the generalized eigenvector corresponding to the least generalized eigenvalue.


The HOSVD-ESPRIT algorithm based on tensor model can be adopted, when the array has translation invariant property; i.e., the array consists of multiple spatial matching subarrays. This section refers to the idea in [10] to derive the HOSVD-ESPRIT model based on the vector-sensor arrays. Assume that the array is obtained by linear translations of the reference subarray, which is obtained by selecting a single vector sensor. Construct the spatial steering vector as follows:where represents the th spatial translation vector, is the number of the th spatial subarray, , and . Therefore, the steering tensor can be constructed as follows:

Furthermore, the tensor-output of snapshots can be denoted aswhere . According to (22), the HOSVD expression of the tensor can be obtained:where the core tensor of same size as . Use to represent the truncated tensor of [10]. The low-rank approximation of can be expressed as

Let and ; denotes the second row of matrix , where . Additionally, if , there is . Furthermore, the parameters estimation of the th mode can be expressed aswhere , , , , and represents the selection matrix for th mode. If the array is uniformly spaced in the th mode, the maximum overlap case can be constructed as follows:

C. ALS-Based CPD for the Third-Order Tensor

The approximate factor matrices of the third-order tensor can be obtained by a solution of the optimization problem:

According to (3), optimization problem (C.1) can be rewrite as a quadratic form:

Furthermore, we can obtain the factor matrix by the following approximate form:

In a similar way, the update methods of the factor matrices and can be derived as

In order to obtain the optimal solution of the problem (C.1), it is necessary to update the matrices , , and through multiple iterations until (C.1) converges to the desired state.

D. Cramér-Rao Lower Bound for the Vector-Sensor Array

Consider the situation described by (14) and establish its matrix form:where and . The unknown parameters vector can be denoted aswhere represents the unknown parameters vector of the th source, . Assume that matrix in (D.1) is column full rank matrix, and the Jacobian is also full rank. Further, setwhere denotes the number of elements in vector . The main concern of this paper is to investigate the performance of estimating in (D.1) from .

To simplify the after-mentioned expression of the Cramér-Rao lower bound, two intermediate matrices are constructed as follows:where denotes the covariance of the signal matrix ,   is the noise variance, and symbolizes a unit matrix. The Cramér-Rao lower bound of the unbiased estimation of iswhere symbolizes a matrix with all entries equal to one, . Assume that with dimension is the -th block entry of the matrix ; the block trace operator , the block transpose operator , the block Kronecker product , and the block Hadamard product are defined as follows.

Definition D.1. Block trace operator is

Definition D.2. Block transpose is

Definition D.3. Block Kronecker product is

Definition D.4. Block Hadamard product is

Data Availability

The simulation data used to support the findings of this study are available from the corresponding author upon request.

Conflicts of Interest

We declare that there are no conflicts of interest regarding the publication of this paper.


This work is supported by The National Key Research and Development Program of China under Grant no. 2018YFB0505104. This work is also supported by Shenzhen Science and Technology Innovation Committee of Basic Research Projects under Grant no. JCYJ20170306154016149.