Abstract
This paper addresses the issue of blind multiple-input multiple-output (MIMO) demodulation of communication signals, with time-varying channels and in an interception context. A new adaptive-blind source separation algorithm, which is based on the implementation of the Multimodulus cost function by analytical methods, is proposed. First a batch processing analysis is performed; then an adaptive implementation of the (Analytical Multi-Modulus Algorithm) AMMA and its simplified version named (Analytical Simplified Constant Modulus Algorithm) ASCMA is detailed. These algorithms, named adaptive-AMMA and adaptive-ASCMA, respectively, are compared with the adaptive (Analytical Constant Modulus Algorithm) ACMA and the MMA (Multi-Modulus Algorithm). The adaptive-AMMA and adaptive-ASCMA achieve a lower residual intersymbol interference and bit error rate than those of the adaptive-ACMA and MMA.
1. Introduction
The context of this study is the Multiple-Input Multiple-Output (MIMO) interception, and more specifically blind demodulation in reception. The channel is supposed to be not stationary, so the methods used to demodulate the received signals must be fast. Such systems exist, among them are the Blind Source Separation (BSS) implemented with analytical methods.
BSS algorithms require few information about emission signals in order to process and recover the transmitted symbols (sources). This type of algorithms estimate an equalization matrix from received signals to obtain the transmitted sources. To find this matrix, different cost function can be minimized, like the Multimodulus (MM) [1, 2], the Constant Modulus (CM) [3, 4], the Simplified Constant Modulus (SCM) [5], and the Multiuser Kurtosis (MUK) [6]. The algorithms in [1–6] use the stochastic gradient to minimize these cost functions, so they are slow to converge and with time-varying environments they are less accurate. So, analytical methods such as the adaptive Analytical Constant Modulus Algorithm (ACMA) [7–10], which converges quickly, should be used with time-varying channels.
Another particularity of the BSS algorithms is their ability to restore, in the output, the sources with a certain permutation and a certain rotation or phase. Moreover, when a frequency offset is present, the constellations in output of adaptive-ACMA and ACMA spin. In order to avoid that and track the frequency offset jointly to the blind source separation, an MM or SCM cost functions can be used. To our knowledge, the MM and the SCM cost functions associated with analytical methods have never been tried so far (we have presented these works partly in conference [11]). So, in this paper an adaptive Analytical Multimodulus Algorithm and its simplified version named the adaptive Analytical Simplified Constant Modulus Algorithm are proposed.
The performances of the adaptive-AMMA and adaptive-ASCMA are compared with those of the adaptive-ACMA and MMA. With time-selective channels, the adaptive-AMMA and the adaptive-ASCMA allow to obtain a higher signal-to-interference noise ratio (SINR) and a lower bit error rate (BER), than those obtained with the adaptive-ACMA and the MMA. Moreover, both proposed methods perform a source separation and a carrier phase recovery, unlike the adaptive-ACMA.
First of all, in this paper, the system model used here is described. Then a batch analysis of the blind source separation algorithm AMMA is performed, and an adaptive form of this algorithm and of its simplified version is presented. Finally, simulation results of the proposed algorithms are compared to those obtained with the adaptive-ACMA and MMA in time-varying envrionments.
The Kronecker and Khatri-Rao products have the following properties:
,,,.2. System Model
A MIMO system with transmitters and receivers is considered. In this system, the information symbols are transmitted through the time-varying MIMO channel at the time . Taking into consideration a carrier frequency offset , the vector of the received signals can be expressed as
where
and
(i) is an representing the baseband received symbols without any time oversampling,(ii) is an vector representing the transmitted signals (sources). The sources are independent, identically distributed (i.i.d.), and mutually independent with zero mean and unit variance,(iii) is an vector representing the additive Gaussian noise,(iv) is an unknown, full column rank, instantaneous linear and time-selective MIMO channel,(v) is the baud rate of the transmitter.3. The AMMA: A BSS Analytical Method so as to Minimize the MM Cost Function
3.1. Hypotheses
The BSS can be performed under the following hypotheses.
(i).(ii) has i.i.d. complex components, is unitary unlike the received signals , and must be prewhitened before applying the BSS.(iii)The noise is additive, white, and Gaussian with a zero mean and a covariance matrix .(iv)The sources are zero mean discrete-time sequences, with a covariance matrix , and they must be mutually independent at a given time and identically distributed.3.2. BSS's Principle
The BSS objective is to find a matrix such as
where the output vector sequence contains accurate estimates of the source signals . Each column of allows to find one source. So, all the columns of must be orthogonal with one another, in order to obtain independent sources. The matrix channel must be unitary, otherwise the received signals will have to be prewhitened. Generally is not unitary, so we will consider a prewithening operation which also reduces the dimension of from to . An underscore is used here to denote pre-whitened variables. Thus, let be the prewhitened received signals, where is an prefilter. Like in [8], the prewhitening algorithm in [12] is used. The initial problem to find is then replaced by finding a matrix :
where and are a global separation matrix. Unfortunately, the calculated matrix separates the sources, except for a possible permutation of and an arbitrary phase on each source signal. The different BSS cost functions, described in Section 3.3, are indeed minimized with any permutation and any rotation on sources. Thus after separation, separator outputs present the following form:
where is a diagonal matrix and is a permutation matrix which represent the arbitrary phase and the permutation, respectively.
In order to find , several cost functions can be used. The next section recalls the well-known (Constant Modulus) CM cost function and (Multimodulus) MM cost function.
3.3. Cost Function
3.3.1. CM Cost Function
Many telecommunication signals have a constant modulus (PSK, 4-QAM), normalized to one generally. If the samples of one source are shown on a complex plan, every samples are on an unitary circle, but it is not true if we show a linear combination of several sources. It is this property that cost function CM exploits. Thus, the algorithms using the CM cost function make use of the sources' constant-modulus property:
where represents the th row of and . This criterion forces BSS outputs to be on a circle with radius .
3.3.2. MM Cost Function
The cost function is composed of two cost functions; one cost function for real part and one for imaginary part of the equalizer output :
The error estimation for the real and imaginary parts is done jointly, and so the cost function uses implicitly the phase of the equalizer output. Thus, the carrier phase recovery is also accomplished with blind equalization, while the equalizer output signal constellation, obtained when using the CM cost function, suffers from an arbitrary phase rotation. So, the MM cost function is most interesting and used later in this paper.
3.3.3. SCM Cost Function
The SCM cost function [5] consists in projection of the equalizer outputs on one dimension (either real or imaginary part). This cost function presents a lower computational complexity compared to the CM and MM cost functions:
3.4. Implementation of Cost Function
Several methods exist to implement cost functions. Among them are the stochastic gradient descent and analytical methods.
The algorithms using stochastic gradient are adaptive algorithms which are slow to converge. Furthermore, in varying environments, subsequent signals tend to be less accurate because of a maladjustment. These algorithms, using a stochastic gradient in order to minimize the cost functions CM and MM, are named, respectively, (Constant Modulus Algorithm) CMA [3, 4] and (MultiModulus Algorithm) MMA [1].
Analytical methods exist in both the adaptive form and batch form. With the adaptive form, the algorithm converges quickly, so it is an interesting option in varying environments.
Our choice is to use, for varying environments, an adaptive analytical method so as to minimize the MM cost function. So, a new algorithm named AMMA and which uses this association, is proposed, once a batch analysis has been performed.
3.5. The Analytical Multimodulus Algorithm
The following cost function should be minimized:
In order to bring out the real and imaginary parts , this last term is written as
where with .
The matrix has the following properties:
Using Kronecker product properties and taking into account [9], and can be written as
3.5.1. The Batch Analysis
Let us denote a fixed number of received signals, the rows of and are then, respectively, stacked in the matrix and .
Using the Khatri-Rao product (), and can be written as and . Let us define the vectors and , the cost function is then given by
with .
The linear system can be rewritten by using a QR decomposition [7, 8].
We have , an unitary matrix, and such as .
After doing the QR factorization of and and applying to and , we obtain
with the column vector composed of ones, a column vector composed of zeros, and two rows vectors and and two matrices, then
The cost function becomes
where and allow to avoid the trivial solution . By squaring (14), the expressions of vectors and , and of matrices and are
By using the definition of the estimated covariance matrices of and , named, respectively, and , and by using the property of the Kronecker product, we obtain the following equalities:
Thus
The criterion is minimized when , and . Afterwards, we consider , as constraints. Since received signals are pre-whitened and thus
Moreover, . Thus, , and to a scaling which is not important the optimization problem becomes
In the appendix, it is shown that vectors , which minimize the function under the hypotheses and , are in the image of (c.f. Theorem 2 and Lemma 2). Moreover, the minimization of is equivalent to minimize (c.f. Theorem 1 and Lemma 1).
Now, we will see how to create vectors that present a Kronecker structure, in the image of .
3.5.2. The Batch Method
The vectors which minimize are eigenvectors associated with the smallest eigenvalues. However, these vectors have not necessarily a Kronecker structure. Thus, in order to obtain the vectors in the image and to have a Kronecker struture, we search a full rank matrix that relates the two bases such so that
with and the and matrix, respectively. Using the property of the Kronecker product, this minimization problem can be rewritten as
with
(i) the th vector of the matrix,(ii) is a diagonal matrix,(iii).Please note that if the matrix was square, the minimization of (23) would be equivalent to a joint diagonalization. In order to obtain the matrix, let us suppose in first phase that it is indeed square. Thus, with the aid of joint diagonalization of matrices, a matrix is obtained. This diagonalization can be done by using the algorithm described in [13]. The columns of the matrix are dependent; the matrix is extracted by keeping the columns of linearly independent.
Now, this algorithm is implemented thanks to an adaptive method. So, we must estimate in an adaptive manner, , , and the eigenvectors of associated with the smallest eigenvalues. The methods used are described in the next section. Then the adaptive tracking of the minor subspace of is performed thanks to the use of the NOOJA algorithm. Finally the adaptive update of the joint diagonalization is realized.
3.5.3. The Adaptive Method
The function will be minimized, at each time , but iteratively
The matrices and were defined by a batch of N samples. In this adaptive method they are converted into an exponential window (-scaling). The details are shown in [8].
The following method used to estimate iteratively the matrix and is a generalization of the Van der Veen method [9]
Then, and are, respectively, the estimations of the autocorrelation matrices for the signals and , with and
where
We saw that the solution of this minimization is in the image of . In order to obtain these vectors in an adaptive manner, a subspace tracking algorithm is used, associated with an adaptive joint diagonalization.
Subspace Tracking
Used alone, the (Normalized Orthogonal Oja) NOOJA algorithm of [14] is utilized to track minor subspace spanned by the eigenvectors relating to the smallest eigenvalues. But, when associated with an adaptive joint diagonalization, as described in the next section, the NOOJA algorithm tracks the minor subspace spanned by the eigenvectors in the image of relating to the smallest eigenvalues which are different from 0. This algorithm extracts adaptively the minor subspace relating to the input signal 's autocorrelation matrix , by maximization the following cost function:
In our problem, the minor subspace of is searched. So, in order to use the NOOJA algorithm, must be expressed as . But since the function (cf. (25) and (26)) is not bilinear, the vector is unknown. Instead of searching the minor subspace of , the minor subspace of both and will be searched. Since can not be expressed as
The matrices and are similar. They have therefore the same image. In order to find the minor subspace spanned by the eigenvectors and relating to the smallest eigenvalues different from 0 for both and , we propose to maximize the following cost function:
Since and have the same image, the first term of the cost function is sufficient to obtain convergence with a constant channel. Using only the matrix is equivalent to minimize the cost function Simplified CM [5] in analytical way. This simplification reduces the complexity of the adaptive AMMA and we name it the adaptive Analytical Simplified CMA (ASCMA). Please note that the ASCMA batch analysis is the same with the AMMA batch analysis performed in Section 3.5.1. However, the study is solely about rather than .
The maximization of may be achieved by using the gradient-descent technique:
where
The time instant has been removed to simplify the equations.
By using , the cost function can be simplified as
By replacing (32) in (34) and by doing the gradient of with respect to , the optimal stepsize is obtained
Keeping in mind that is orthogonal at each iteration, that is, , and replacing and by the estimate and , respectively, we can write
where , and . The optimal step size becomes
The instant has been removed in the 's expression in order to simplify the equation.
The updating of is performed as
where the matrix is obtained thanks to a Gram-Schmidt orthogonalization of . The algorithm obtained is described in Algorithm 1.
Once the eigenvectors relating to the smallest eigenvalues of both and are obtained, the constraint must be satisfied. In other words, must have a Kronecker structure and must be in the image of .
| ||||||||||||||||||||||||||||
Adaptive Update of the Joint Diagonalization
The method used to satisfy the constraint is done as in [9]. is computed and regarded as the current estimate of the subspace basis (). Using this basis, the subspace update is performed by using the “Normalized Orthogonal OJA” giving . The last step of this algorithm is the mapping of the columns of to a Kronecker-product:
Then a power iteration [15] is applied. This iteration can be written as
At the next iteration, the matrix in the image, which has a Kronecker structure, is obtained with . Then, is transmitted at the NOOJA algorithm, allowing to vectors obtained in output NOOJA to be in the nullspace and to have a Kronecker struture.
Estimated Symbols
In order to obtain the transmitted symbols, the complex matrix is rebuilt from
In order to prevent the extraction of the same source many times, is reconditioned after each update. The technique used here is the computation of a singular value decomposition of : . Then the singular values of , which are smaller than 0.5, are replaced by 1
Then, the matrix is applied on the pre-whitened received signals:
The suggested adaptive-AMMA algorithm is described in Algorithm 2.
4. Simulation Results
The channel used in simulations is a time-selective channel; that is, it is time varying but flat in frequency. The transmitted symbols, resulting from a 4-QAM constellation, are sent through a MIMO Rayleigh fading channel. Two transmitters and four receivers are used ( and ). The performances of proposed adaptive-AMMA are compared with those of the adaptive-ASCMA, the adaptive-ACMA and the MMA, in terms of Bit Error Rate (BER) and Signal to Interference plus Noise Ratio (SINR).
In the simulation results of the Figure 1, we considered the case of a carrier frequency offset () with dB. The equalizer output constellation of the adaptive-ACMA (Figure 1(e)) spins because of the carrier frequency offset. Normally, the MMA (Figure 1(f)) allows to track the carrier frequency offset, but in this case the carrier phase is too high to obtain a good tracking. Only the proposed adaptive-AMMA and adaptive-ASCMA can remove the frequency offset (Figures 1(c) and 1(d)), without the use of a carrier tracking loop.
(a) Received signals
(b) Transmitted signals
(c) Adaptive-AMMA
(d) Adaptive-ASCMA
(e) Adaptive-ACMA
(f) MMA
Figure 2 shows the average BER as a function of SNR (with the term “average” being here a time average), and this for two sources and 1000 runs, each one presenting different randomly varying channels and different noise sequences. The value of the Doppler shift is . The simulation here was realized without any carrier frequency offset.
The lowest BER on Figure 2 is obtained with the proposed adaptive-AMMA algorithm, followed by the adaptive-ASCMA, then by the adaptive-ACMA and the highest BER is obtained with the MMA.
So the proposed adaptive-AMMA algorithm allows to obtain a 4 dB gain in SNR terms compared with that obtained by the adaptive-ACMA with a value of BER equal to .
Figure 3 illustrates the convergence speed of the adaptive-AMMA. The path gain of the channel and its estimate obtained by the adaptive-AMMA are represented on the same figure. We can notice that estimation of adaptive-AMMA follows the variations of the channel.
From the results of the Figure 4, we can see that the performances of the proposed adaptive-AMMA and the adaptive-ASCMA are very closed in residual SINR terms. Moreover, these two algorithms perform higher residual SINR than the adaptive-ACMA and the MMA.
5. Conclusion
In this paper, we have proposed two new BSS algorithms, named adaptive-AMMA and adaptive-ASCMA, build on the Multimodulus criteria, the Simplified Constant Modulus respectively, and an analytical processing method. The second algorithm is a simplified version of the adaptive-AMMA and it accomplishes performances results closed to the adaptive-AMMA. These algorithms perform better results in BER and SINR terms compared with the adaptive-ACMA and the MMA algorithms. In addition to the blind separation, they can perform a carrier phase recovery simultaneously unlike the adaptive-ACMA.
Appendix
Theorem 1. The matrices and are symmetric and semidefine positive.
Proof. Let us consider a symmetric matrix, so and thus is a projection. Then . Now, we have Similarly Since and , then and are symmetric, semidefine positive and so too. Furthermore, .
Lemma 1. and so the optimization problem becomes
Afterward, let us define the nullspace and the image of as
and are subspaces of .
Property 1. where denotes the direct sum, that is, and and are subspaces of .
Let us assume now, for the sake of simplicity, that and , .
Afterward, the following affirmations will be proven. The vectors which minimize are in the nullspace of , but the sole vectors verify the constraint , yet is a trivial solution to be avoided. However, the vectors are in the image of .
Definition 1. By construction, the sole dependence obtained on the columns of is , where represents the th column of . So, . A basis of the nullspace is defined as
where are vectors in the standard basis of .
The are linearly independent and verify .
Thanks to the following rank-nullity theorem:
and since the dimension of is equal to 6, the dimension of is equal to 10. And is a basis of image:
The are linearly independent; verify and we can verify the Property 1. The and the are linearly independent, then is a basis for .
Now, the following theorem can be proved.
Theorem 2. In the nullspace of , the vector is the sole vector to have a Kronecker structure.
Proof. Let , . Supposing ,
We have(i),(ii), , , , , .So, is a skew-symmetric matrix.
and with are equal to
is a symmetric matrix and its diagonal values are positive, while is a skew-symmetric matrix. So, vectors do not have a Kronecker structure, except .
So, vectors are not in the nullspace of , but vectors are in the image of . Using property 1, we deduce the following lemma.
Lemma 2. Let us consider , the set of vectors in , which have a Kronecker structure: , then .
Notations
The following notations are adopted:
| : | Complex conjugation |
| : | Matrix or vector transpose |
| : | Matrix or vector complex conjugate transpose |
| : | Real part of complex |
| : | Imaginary part of complex |
| : | Vector of all 1s |
| : | An matrix with zeros components |
| : | The identity matrix |
| : | Expectation operator |
| : | A diagonal matrix constructed from a vector |
| : | Stacking of the columns of into a vector |
| : | Kronecker product |
| : | Khatri-Rao product: |
| : | The Frobenius norm |