一种多载波CDMA信号扩频码及信息序列盲估计方法

摘要

本发明公开了一种多载波CDMA信号扩频码及信息序列盲估计简化方法，通过建立频域扩频的信号模型，考虑时变的信道环境，结合分段建模的思想，采用序贯蒙特卡罗的抽样粒子来近似建立的联合后验概率分布，并根据抽样值和重要性权值大小进行参量状态估计；同时合理的设计执行流程，修正原有的迭代计算重要性权值步骤，并充分利用空时分组码和多载波调制的正交属性，简化方法的实现过程。本发明较好地满足了时变信道环境下具有空时分组码的多载波CDMA信号扩频码及信息序列的快速提取要求，大大提高了扩频参数的估计精度，并且具有较低的计算复杂度。本发明可直接应用于非协作扩频通信系统，也可用于相应的软件无线电等系统。

说明书

技术领域

本发明属于信号处理领域中非协作通信信号处理技术，具体是指一种具有空时分组码的多载波CDMA信号扩频码及信息序列盲估计方法。

背景技术

随着DS-CDMA技术的不断发展，传统的扩频技术因为存在处理增益和数据率之间的矛盾，从而限制了其在高速数据传输方面的应用。为了更好的满足卫星通信的需求，国外现有的卫星广泛采用WCDMA通信体制，并在扩频通信中引入了多载波(MC)直接序列扩频技术，同时采用空时分组码(STBC)进行编码，以提高数据的传输速率和增加系统的容量。因此，在非合作通信条件下，对具有空时分组码的多载波CDMA信号进行扩频码及信息序列等参数估计具有重要的军事意义。

目前，在非合作信号处理领域，对具有空时分组码的多载波CDMA信号的处理研究仍处于起步阶段。在《IEEE Journal on Selected Areas in Communiations》杂志2001年19期“Bayesian Monte Carlo MultiuserReceiver for Space-Time Coded Multi-carrier CDMA Systems”一文中，Yang Zigang在合作条件下采用吉布斯抽样的方法对该系统信息序列和信道参数进行了有效的估计，得到了较好的结果；但针对非合作情况，现有的研究成果主要是针对多载波CDMA进行的，不能直接植入应用于STBC-MC-CDMA系统中，且这些方法对低信噪比的适应能力和计算复杂度有待进一步改善。2010年Bangwon Seo采用最小均方误差自适应滤波的方法进行处理，其收敛速度较快，但对信噪比的适应能力不强。在此基础上，2012年Tsui-Tsai Lin在(IEEE Transactions on Vehicular Technology)上公开了一种针对STBC-MC-CDMA系统的相关矩阵方法，该方法能很好的适应低信噪比的需求，且计算量相对较低，其不足之处是：收敛速率较慢，鲁棒性较差。

由此可以看出，已有的方法还不能满足STBC-MC-CDMA系统的需要，同时考虑时变的信道环境，还需研究一种更有效的扩频码及信息序列盲估计简化方法。

发明内容

本发明所要解决的技术问题是，针对现有技术的不足，提出了一种基于序贯蒙特卡罗(SMC)分析的STBC-MC-CDMA信号扩频码及信息序列盲估计方法，同时考虑时变信道环境等非线性因素的影响，它可以较好地满足非协作扩频通信中扩频码及信息序列的快速提取要求，大大提高了扩频参数的估计精度，并且具有较低的计算复杂度。本发明可直接应用于非协作扩频通信系统，也可用于相应的软件无线电等系统。

为解决上述技术问题，本发明是通过以下技术方案实现的：通过建立频域扩频的信号模型，考虑时变的信道环境，结合分段建模的思想，采用SMC的抽样粒子来近似建立的联合后验概率分布，并根据抽样值和重要性权值大小进行参量状态估计；同时合理的设计执行流程，修正原有的迭代计算重要性权值步骤，并充分利用空时分组码和多载波调制的正交属性，简化方法的实现过程。

所述的扩频码及信息序列盲估计方法主要用于STBC-MC-CDMA信号的扩频码及信息序列盲估计，且信号模型的建立考虑时变的信道环境，假设第k个用户的第m帧数据在第j个天线上经过多径信道传播，其表达式为

$>h_{k,m,j}\left(t\right)=\underset{l=0}{\overset{L-1}{\Sigma }}{\beta }_{k,m,j,l}\delta (t-\frac{l}{\mathrm{\Delta l}})---\left(1\right)$>

式中：δ(·)为狄拉克函数，L为传播的路径数，β_k,m,j,l为第k个用户的第m帧数据在第j个天线上第l条路径的复衰落幅度，Δ_f为多载波系统的带宽。在经过离散傅里叶变换之后，信道在第n个子载波上的频域响应为

$>H_{k,m,j}\left(n\right)\stackrel{\Delta }{=}{H}_{k,m,j}\left(n{\Delta }_f\right)=\underset{l=0}{\overset{L-1}{\Sigma }}{\beta }_{k,m,j,l}\exp (-\mathrm{j\pi nl}/N)=w_f^H\left(n\right)h_{k,m,j}---\left(2\right)$>

式中：$>h_{k,m,j}\stackrel{\Delta }{=}{[{\beta }_{k,m,j,0},{\beta }_{k,m,j,1},...,{\beta }_{k,m,j,L-1}]}^T$>为多径的时域响应，令$>h_{m,j}\stackrel{\Delta }{=}{[h_{1,m,j}^T,h_{2,m,j}^T,...,h_{K,m,j}^T]}^T,$>为相应的离散傅里叶变化系数，即

$>w_f\left(n\right)\stackrel{\Delta }{=}{[1,\exp (-j2\mathrm{\pi n}/N),...,\exp (-j2\mathrm{\pi n}(L-1)/N)]}^H---\left(3\right)$>

假设信道在一帧的STBC时间内是非时变的，而在不同的M内，将时变的衰落幅度建模为一阶AR模型，即

h_m,j＝Fh_m-1,j+w_m              (4)

式中：F＝αI_KL,w_m为零均值，协方差矩阵为的复高斯分布。

于是在接收端取M帧的STBC数据进行处理(此时包含2M帧的时域数据，因为在一个分组内传送2帧的数据)，可得用户总数为K的离散接收信号为

式中：T为交错函数，在一般的系统中，T为(N×N)维的初等矩阵，v_m为零均值协方差矩阵为的复高斯噪声。且有

$>C_{k,j}\stackrel{\Delta }{=}\mathrm{diag}{\left\{c_{k,j}\right\}}_{N\times N},c_{k,j}={[c_{k,j}{\left(1\right)}^T,...,c_{k,j}{\left(P\right)}^T]}^T,j=\mathrm{1,2}---\left(6\right)$>

$>B_{k,m,j}\stackrel{\Delta }{=}\mathrm{diag}\{b_{k,m,j}\otimes I_G\}=\mathrm{diag}{\{b_{k,m,j}\left(1\right)I_G,...,b_{k,m,j}\left(P\right)I_G\}}_{N\times N},j=\mathrm{1,2}---\left(7\right)$>

$>W_f\stackrel{\Delta }{=}{[w_f\left(0\right),w_f\left(1\right),...,w_f(N-1)]}_{N\times L}^H---\left(8\right)$>

式中：代表Kronecker乘积，b_k,m,j,j＝1,2为第m帧STBC时刻需发送的符号向量，c_k,j为分配给第j个天线的扩频码向量。

所述的分段建模思想指的是：将接收信号按照不同的符号路数和天线数目进行分段处理，分别对p＝1,…,P和j＝1,2取不同值时采用SMC方法进行迭代抽样，其信号模型为

$>\left(\begin{array}{c}{y}_m(p,j)=\underset{k=1}{\overset{K}{\Sigma }}{X}_{k,m}(p,j)W(p,j)h_{k,m,j}+v_m(p,j)\\ m=1,...,M.P=1,...,P.j=\mathrm{1,2}\end{array}\right)---\left(9\right)$>

式中：X_k,m(p,j)为从矩阵X_k,m中抽取的(G×G)维子矩阵，该子矩阵由包含元素b_k,m,j(p)的所有行和列组成；同理y_m(p,j)、..、h_k,m,j和v_m(p,j)也为从对应矩阵中抽取组成的子矩阵。由于对所有的用户K，交错函数T是相同的，于是将上式写成矩阵相乘的形式可得

y_m(p，j)＝X_m(p，j)Ψ(p，j)h_m，j+v_m(p，j)

                                       (10)

m＝1，…，M.p＝1，…，P.j＝1，2式中：h_m,j为KL维列向量，且有

$>X_m(p,j)\stackrel{\Delta }{=}{[X_{1,m}(p,j),...,X_{K,m}(p,j)]}_{G\times \mathrm{KG}}---\left(11\right)$>

为表述方便，作如下定义

b_m(p,j)＝[b_1,m,j(p),b_2,m,j(p),…,b_K,m,j(p)]           (13)

B_m(p,j)＝[b₁(p,j),b₂(p,j)…,b_m(p,j)]              (14)

c(p,j)＝[c_1,j(p)^T,c_2,j(p)^T,…,c_K,j(p)^T]^T           (15)

Y_m(p,j)＝[y₁(p,j),y₂(p,j),…,y_m(p,j)]           (16)

$>H_{m,j}\stackrel{\Delta }{=}{[h_{1,j}^T,h_{2,j}^T,...,h_{m,j}^T]}^T---\left(17\right)$>

本发明的研究目标是：在p,j分别取不同值时，根据SMC方法，从联合后验分布p(B_m(p,j),c(p,j),H_m,j|Y_m(p,j))中序贯抽样得到粒子1≤m≤M，从而估计出对应的扩频码和信息序列向量。

根据贝叶斯准则，将未知参量的后验分布函数分解如下：

p(B_m(p，j)，c(p，j)，H_m，j|Y_m(p，j))

                                         (18)

＝p(H_m，j|B_m(p，j)，c(p，j)，Y_m(p，j))p(B_m(p，j)，c(p，j)|Y_m(p，j))

从上式可以看出，对未知参量的估计可以分成两部分进行，对于等号右边第一项可以采用卡尔曼滤波方法对H_m,j的均值和协方差进行迭代更新，第二项为低维的后验分布概率，采用SMC方法进行序贯抽样，令Z_m(p,j)＝(B_m(p,j),c(p,j))，z_m(p,j)＝(b_m(p,j),c(p,j))。

所述的简化实现过程指的是：结合卡尔曼滤波方法的执行流程，充分利用空时分组码和多载波调制的正交属性，简化迭代步骤，大大减少方法的计算复杂度，具体如下

由式(4)和(10)可知，在第m帧STBC时间内，h_m,j服从复高斯分布，即

$>p\left(h_{m,j}\right|Z_m^{\left(i\right)}(p,j),Y_m(p,j))~N_c(h_{m|m,j}^{\left(i\right)},P_{m|m,j}^{\left(i\right)})---\left(19\right)$>

式中的均值和方差可以通过卡尔曼滤波方法，在给出和时进行迭代更新，其过程为：

$>h_{m|m-1,j}^{\left(i\right)}={\mathrm{Fh}}_{m-1|m-1,j}^{\left(i\right)}---\left(20\right)$>

$>P_{m|m-1,j}^{\left(i\right)}={\mathrm{FP}}_{m-1|m-1,j}^{\left(i\right)}{F}^H+{\Sigma }_w---\left(21\right)$>

$>{\hat{y}}_m^{\left(i\right)}(p,j)=X_m^{\left(i\right)}(p,j)\Psi (p,j)h_{m|m-1,j}^{\left(i\right)}$>

$>{\upsilon }_{m,j}^{\left(i\right)}=y_m^{\left(i\right)}(p,j)-{\hat{y}}_m^{\left(i\right)}(p,j)---\left(23\right)$>

$>P_{\upsilon ,m,j}^{\left(i\right)}={\Sigma }_v+\left[X_m^{\left(i\right)}(p,j)\Psi (p,j)\right]P_{m|m-1,j}^{\left(i\right)}{\left[X_m^{\left(i\right)}(p,j)\Psi (p,j)\right]}^H---\left(24\right)$>

$>G_{m,j}^{\left(i\right)}=P_{m|m-1,j}^{\left(i\right)}{\left[X_m^{\left(i\right)}(p,j)\Psi (p,j)\right]}^H{\left[P_{\upsilon ,m,j}^{\left(i\right)}\right]}^{-1}---\left(25\right)$>

$>h_{m|m,j}^{\left(i\right)}=h_{m|m-1,j}^{\left(i\right)}+G_{m,j}^{\left(i\right)}{\upsilon }_{m,j}^{\left(i\right)}---\left(26\right)$>

$>P_{m|m,j}^{\left(i\right)}=P_{m|m-1,j}^{\left(i\right)}-G_{m,j}^{\left(i\right)}{P}_{\upsilon ,m,j}^{\left(i\right)}{\left(G_{m,j}^{\left(i\right)}\right)}^H---\left(27\right)$>

式中：为状态预测值，为状态预测值的协方差矩阵，为预测的观测值，为新息过程，为新息矩阵，为卡尔曼增益。

由式(9)和(10)的定义可知

$>\left(\begin{array}{c}{\left[X_m(p,j)\Psi (p,j)\right]}^H\left[X_m(p,j)\Psi (p,j)\right]={\Psi }^H(p,j)X_m^H(p,j)X_m(p,j)\Psi (p,j)\\ =\mathrm{KG}·I_{\mathrm{KL}}\end{array}\right)---\left(28\right)$>

式中用到了STBC的正交属性以及多载波调制的正交属性Ψ^H(p,j)Ψ(p,j)＝KG·I_KL。假设为对角矩阵，即

$>P_{m-1|m-1,j}^{\left(i\right)}={\delta }_{m-1}^{\left(i\right)}{I}_{\mathrm{KL}}---\left(29\right)$>

将上式代入式(21)可得也为对角矩阵，即

$>P_{m|m-1,j}^{\left(i\right)}={\delta }_{m-1}^{\left(i\right)}{\mathrm{FF}}^H+{\Sigma }_w=({\left|\alpha \right|}^2{\delta }_{m-1}^{\left(i\right)}+{\sigma }_w^2)I_{\mathrm{KL}}---\left(30\right)$>

令并根据矩阵的求逆定理，将式(24)新息矩阵的逆写成

$>\left(\begin{array}{c}{\left[P_{\upsilon ,m,j}^{\left(i\right)}\right]}^{-1}={\Sigma }_v^{-1}-{\Sigma }_v^{-1}\left[X_m^{\left(i\right)}(p,j)\Psi (p,j)\right]\\ \times {\left({\left[X_m^{\left(i\right)}(p,j)\Psi (p,j)\right]}^H{\Sigma }_v^{-1}\right[X_m^{\left(i\right)}(p,j)\Psi (p,j)]+{\left(P_{m|m-1,j}^{\left(i\right)}\right)}^{-1})}^{-1}{\left[X_m^{\left(i\right)}(p,j)\Psi (p,j)\right]}^H{\Sigma }_v^{-1}\end{array}\right)---\left(31\right)$>

代入式(28)和(30)可得

$>\left(\begin{array}{c}{\left[P_{\upsilon ,m,j}^{\left(i\right)}\right]}^{-1}=\frac{1}{{\sigma }_v^2}{I}_G-\frac{1}{{\sigma }_v^4}\left[X_m^{\left(i\right)}(p,j)\Psi (p,j)\right]\\ \times {\left(\right[\frac{\mathrm{KG}}{{\sigma }_v^2}+\frac{1}{{\beta }_m^{\left(i\right)}}\left]I_{\mathrm{KL}}\right)}^{-1}{\left[X_m^{\left(i\right)}(p,j)\Psi (p,j)\right]}^H\\ =\frac{1}{{\sigma }_v^2}{I}_G-\left(\frac{{\beta }_m^{\left(i\right)}}{{\sigma }_v^4+\mathrm{KG}{\sigma }_v^2{\beta }_m^{\left(i\right)}}\right)\left[X_m^{\left(i\right)}(p,j)\Psi (p,j)\right]{\left[X_m^{\left(i\right)}(p,j)\Psi (p,j)\right]}^H\end{array}\right)---\left(32\right)$>

于是结合式(25)和(27)可得为

$>\left(\begin{array}{c}{P}_{m|m,j}^{\left(i\right)}=P_{m|m-1,j}^{\left(i\right)}-P_{m|m-1,j}^{\left(i\right)}{\left[X_m^{\left(i\right)}(p,j)\Psi (p,j)\right]}^H{\left({\left[P_{\upsilon ,m,j}^{\left(i\right)}\right]}^{-1}\right)}^H\left[X_m^{\left(i\right)}(p,j)\Psi (p,j)\right]{\left(P_{m|m-1,j}^{\left(i\right)}\right)}^H\\ =\left(\frac{{\beta }_m^{\left(i\right)}{\sigma }_v^2}{\mathrm{KG}{\beta }_m^{\left(i\right)}+{\sigma }_v^2}\right)I_{\mathrm{KL}}={\delta }_m^{\left(i\right)}{I}_{\mathrm{KL}}\end{array}\right)---\left(33\right)$>

由式(30)和(33)可知，在取为对角矩阵时，和也为对角矩阵，于是式(20)～(27)可以简化为

$>h_{m|m-1,j}^{\left(i\right)}={\mathrm{Fh}}_{m-1|m-1,j}^{\left(i\right)}---\left(34\right)$>

$>{\beta }_m^{\left(i\right)}={\left|\alpha \right|}^2{\delta }_{m-1}^{\left(i\right)}+{\sigma }_w^2---\left(35\right)$>

$>{\mu }_m^{\left(i\right)}=\frac{{\beta }_m^{\left(i\right)}}{\mathrm{KG}{\beta }_m^{\left(i\right)}+{\sigma }_v^2}---\left(36\right)$>

$>h_{m|m,j}^{\left(i\right)}=(1-{\mu }_m^{\left(i\right)}\mathrm{KG})h_{m|m-1,j}^{\left(i\right)}+{\mu }_m^{\left(i\right)}{\left[X_m^{\left(i\right)}(p,j)\Psi (p,j)\right]}^H{y}_m^{\left(i\right)}(p,j)---\left(37\right)$>

$>{\delta }_m^{\left(i\right)}=\frac{{\beta }_m^{\left(i\right)}{\sigma }_v^2}{\mathrm{KG}{\beta }_m^{\left(i\right)}+{\sigma }_v^2}---\left(38\right)$>

采用式(34)～(38)进行计算时，在初始时刻取为对角矩阵，在迭代过程中就可以避免计算式(25)中的逆，大大减少了算法的计算量。

所述的修正方法指的是首先对所有的进行迭代处理，估计出各用户的扩频码，再对观测数据进行cholesky分解处理，分别得到各用户的观测处理数据，最终抽样得到具有最大后验概率的状态参量值。

与现有技术相比，本发明的有益效果是：

(1)能够在较低的信噪比条件下适应具有空时分组码的多载波CDMA信号，在不同的扩频长度、不同的粒子数和时变多径信道环境下都具有较高的估计性能，且计算量较低。

(2)能够同时对扩频码及信息序列进行估计，在非协作通信系统中可以省去解扩和解调的过程，极大地简化了系统。

(3)本发明可以完成对中频信号的处理，同时也可以完成对基带信号或解调码流的处理，适应数据类型广泛。

附图说明

图1是本发明所述方法的总体流程图；

图2是MC-CDMA系统的发射结构；

图3是空时分组码的传输结构；

图4是实施例中不同用户时的扩频码估计性能比较；

图5是实施例中不同用户时的信息序列估计性能比较；

图6是实施例中不同扩频长度时的扩频码性能比较；

图7是实施例中不同粒子数时的信息序列估计性能比较；

图8是实施例中时变衰落信道幅度的跟踪抽样均值。

具体实施方式

下面结合附图和实施例对本发明作进一步详细描述。

图1是本发明的总体流程图，本实施例所述方法实现过程如下：

(1)接收中频采样数据或基带数据，估计扩频信号的扩频周期、码片速率和用户个数，并将其传送到扩频码及信息序列估计模块。

(2)初始化数据。对i＝1,…,I，根据B_m(p,j)和c(p,j)所服从的离散均匀先验分布初始化和取$>{\delta }_0^{\left(i\right)}={ϵ}_i,0<{ϵ}_i<1,P_{0|0,j}^{\left(i\right)}={\delta }_0^{\left(i\right)}{I}_{\mathrm{KL}},$>服从均值为0，协方差矩阵为的复高斯分布；令m＝1，$>w_0^{\left(i\right)}(p,j)=1.$>

(3)对于固定的首先采用序贯蒙特卡罗的方法对所有的进行迭代处理，得到c_m(p,j)的估计，并根据序贯抽样的思想计算其次采用cholesky分解方法对观测数据进行处理，分别得到各用户的观测处理数据；然后根据式(34)～(38)对均值和协方差矩阵进行更新；最后抽样得到具有最大后验概率的b_m(p,j)。

(4)计算重要性权值

(5)重抽样。计算有效抽样尺度其中为重要性权值的方差，当I_m小于预先设定的门限值时，进行重抽样，得到新粒子其重要性权值均为$>{\{{\hat{w}}_m^{\left(i\right)}(p,j)=1\}}_{i=1}^I.$>

(6)重复步骤3～步骤5进行多次迭代，最终得到各状态参量估计结果，进而对扩频码和信息序列进行拼接和重构，同时若m＝200则退出，否则m＝m+1，进入下一时刻的计算。

实施例中：用户数K分别为4、6、8、10，子载波数N＝128，每个用户在第p路数据上(p＝1,2)的扩频码为长度G＝64的随机序列，各用户的多径路数均为L＝3，将时变信道衰落幅度建模成一阶AR模型，即h_m,j＝Fh_m-1,j+w_m，取F＝0.999I_KL,w_m为零均值，协方差矩阵为Σ_w＝0.01I_KL的复高斯分布，采用M＝200帧的STBC数据进行处理，迭代过程中取粒子数I＝100。

图4和图5分别为用户数目K为4、6、8、10时，本发明得出的扩频码和信息序列估计性能随输入信噪比变化的曲线。在多用户情况下，输出信噪比和误码率的计算为对所有用户的扩频码和信息序列粒子抽样值求平均得到，且图4中的扩频码输出信噪比是对两个天线的结果求平均后获得。实施例结果表明，在不同用户数目时，本发明能在较低信噪比条件下盲估计出扩频码和信息序列，当信噪比大于-5dB，用户数不超过10时，恢复出的扩频码平均输出信噪比大于13dB，所有用户的信息序列估计误码率均低于10^-2。

图6给出在不同扩频长度条件下，扩频码的估计性能。从图中可以看出，随着G的增加，序列的估计性能随着提高。其主要原因是在迭代过程中，分组的路数越多，信号间产生交叠的概率就越大，从而导致算法的性能降低。

图7给出了G＝64，粒子数I分别为50、100、200时信息序列误码率随输入信噪比变化曲线。由图可以看出，粒子数越大，其估计性能越好，但是计算量也成倍增加，因此在选择粒子数时应折中考虑。

图8给出了当用户数K＝4，SNR＝-6dB，迭代进行到m＝N＝200时，时变衰落信道幅度β_1,m,1,1的跟踪抽样均值。