MIMO系统检测仿真
一、引言
随着无线通信业务的发展,人们对数据率的要求越来越高,而传统通信方式通过使用某些信道编码方法已接近香农极限,要想再提高频谱利用率已经很困难。在这种情况下,多输入多输出(MIMO,Multiple Input Multiple Output)技术由于能同时带来分集增益和空间复用增益,成为未来移动通信系统的有力竞争方案。MIMO通信系统的检测器是MIMO技术实用过程中关键的一个模块,选择一种检测性能好而且便于硬件实现的检测方法是人们追求的目标。
传统的MIMO检查算法主要有:最大似然(ML,Maximum Likelihood)检测算法、迫零(ZF,Zero Forcing)检测算法、最小均方误差(MMSE,Minimum Mean-Square Error)检测算法、V-BLAST(ZF-OSIC)检测算法和基于QR分解的检测算法等。此外,通过把在给定格中寻求最短向量的球形解码思想应用于MIMO系统,形成了MIMO系统的球形解码算法,在保持优良检测性能的同时,大大减小了计算复杂度。本次课程设计主要针对最大似然算法,迫零算法和最小均方误差算法进行仿真和性能仿真比较。
二、MIMO系统
MIMO通信系统可以定义为收发两端分别采用多个天线或阵列天线的无线通信系统。MIMO的多输入多输出是针对多径无线传输信道而言的。
考虑n T根发射天线n R根接收天线的MIMO系统,如下图所示,数据流被分成n T个子数据流,每个子流通过星座点映射后送给发射天线。分别从个发射天线发射出去,再经多径传输信道后由n R个接收天线接收,同时用接收到的信号进行信道估计得到信道参数值,然后通过一定的检测算法处理分解出子信息流。因为n T个发射天线同时发射子信息流,各发射信号只占用同一频带,并未增加带宽,达到提高频谱利用率的目的,同时多个并行空间也实现了更高的数据传输速率。
在接收端的一根天线会收到每根发送天线送出的信号,将所有接收天线收到的符号作为一个矢量12(,,)R
T n x x x x =…,表示,那么x Hs n =+,12(,,)T T n s s s =…,s 是发射信号矢量,H是R T n n ?维的矩阵,其元素,j i h 是发射天线(1,2,,)T i i n =…到接收天线(1,2,,)R j j n =…的信道增益,12(,,)R T n n n n n =…,是各分量独立且都服从),0(2σN 分布
的复白高斯噪声。三、检测算法
(一)ZF 算法
迫零检测是MIMO 系统中常用的检测器,其核心思想是在接收端通过线性变换消除不同天线发射信号间的干扰。将MIMO 系统的信号检测模型改写成如下形式:1122T T n n x h s h s h s n =+++L
其中(1,2,)i T h i n =L 是H 的第i 列。
为了在接收端恢复(1,2,)i T s i n =L 而排除其他分量的干扰,可以使用矢量i w 与x 作内积,其中i w 满足如下条件:
01i j i j w h i j ≠?=? =?
将(1,2,,)i T w i n =L 作为行向量组成一个矩阵ZF W ,显然它应该满足ZF
W H I ?=,所以1()H H ZF W H H H -=(假设H 列满秩),此时发射信号估计值s
$为 1()()H H ZF s H H H Hs n s W n - =+ =+
$ 协方差矩阵为
{}
21()()H H ZF C E s
s s s H H σ-=--= ()$$ 从上面这些式子可以看出,经过迫零检测器之后得到的对发射信号的估计值,完全消除了不同天线发送的数据之间的干扰,在高信噪比条件下有较好的性能。特别地,当噪声项为
0时,严格地有s
s =$。但在低信噪比或者信道矩阵H 接近奇异时,检测性能严重恶化。(二)MMSE 算法
最小均方误差检测则是基于最大化输出信干噪比(SINR ,Signal-Interference -and-Noise Ratio )的考虑,在抑制噪声和消除干扰之间找到一个最佳的平衡点。MMSE 检测的目标是找到估计值s
Wx =$,使其与真实值s 的差异尽可能小。MMSE 的目标函数如下所示
2arg min ||||MMSE W
W E s Wx =- 经过求解得212()H H MMSE W H H I H p σ-=+
,其中2()H E ss p I =,此时估计量的协方差矩阵
为 {}2
21
2()()H H MMSE C E s s s s H H I p σσ-=--= (+)$$
(三)ML 检测算法
最大似然(ML )检测计算接收信号向量和所有可能的后处理向量之间的欧氏距离,并找到一个最小距离,即所有可能的发出的信号矢量中进行搜索:
当所有发射向量等可能性时,ML 算法达到最大后验概率检测的最佳性能,但它的复杂度随调制阶数和发射天线的数量增加而上升,计算度量复杂度随天线数呈指数上升。但因为它具有最佳的性能,尽管计算复杂度比较高,仍将其作为其他检测方法的参考。改进的ML 算法检测可以将ML 度量复杂度降低,但N TX 大于3时复杂度仍然很高。四、实验结果与分析
通过MATLAB仿真,对ML, ZF, MMSE三种算法进行误比特率(BER)性能分析。假设信道是具有丰富散射环境的平坦瑞利信道,发射天线间距和接收天线间距足够大,接收端确知信道状态且能够保持精确同步,各个天线之间的信道参数为零均值单位方差独立分布的复高斯随机变量。采用BPSK调制,天线配置2x2,信噪比范围为0~25dB,发送符号数目为10^6。
对于ML算法,分别考虑 [s1 s2 ] = [+1,-1][+1,+1],[-1,+1],[-1,-1]四种情况,从四组数值中找出最小值并记录所在位置,然后统计错误个数。对ZF算法,首先计算伪逆矩阵 G = inv(H^H*H)*H^H,然后将H^H*H 的矩阵维数记为[nTx x nTx],即[2 x 2],并根据[a b; c d] 的逆1/(ad-bc)[d -b;-c a]求逆,然后根据接收端硬判决解码计算其统计错误。对于MMSE算法我们也是类似地计算W = inv(H^H*H+sigma^2*I)*H^H,其中H^H*H i的维数为[2 x 2],并按照[a b; c d]=1/(ad-bc)[d -b;-c a]求逆,然后根据接收端硬判决解码计算其统计错误。
应用ML, ZF, MMSE三种算法进行信号检测的BER性能曲线仿真图如下:
在同等的情况下,ML检测的性能优于其他两种,MMSE检测的性能次之,ZF检测的性能最差。ZF检测算法会可能带来对高斯噪声的放大,从而影响了检测的准确性,而MMSE检测算法是在ML的基础上试图消除检测算法对噪声的方法,很好的抑制噪声,得到了比ZF检测算法更好的性能。这两种算法均属于线性调制,采用的硬判决的方式。相比于ML检测来说,
ML检测是理论上的最优,在实际的应用中随着天线数目的增加和更多进制调制方式的采用会使得ML检测的计算量成指数增加,设计更复杂。
考虑计算的复杂度问题,性能最优的ML检测算法复杂度是指数形式,算法复杂度随着发射天线数和调制阶数指数增长,算法复杂度为O(S M),线性检测算法ZF和MMSE虽然性能较差,但计算复杂度主要集中在矩阵求逆,复杂度为O(M3)。
五、总结
MIMO系统接收端收到的信号在时间和频率上式重叠的,可能会发生码间干扰,检测难度远远高出传统的单输入单输出的系统,如何在接收端将发射信号分离并正确检测发射信号是MIMO系统的关键问题。这次课程设计比较了几种基本的MIMO检测算法,进行了误比特率性能分析和复杂度分析,认识了主要算法的优缺点。
六、部分代码
发送端
N = 10^6; % 发送的符号数目
Eb_N0_dB = 0:25; % 信噪比范围
nTx = 2;
nRx = 2;
for ii = 1:length(Eb_N0_dB)
% 发送端
ip = rand(1,N)>0.5; % 等概率产生0和1
s = 2*ip-1; % BPSK 调制 0 -> -1; 1 -> 0
sMod = kron(s,ones(nRx,1)); %
sMod = reshape(sMod,[nRx,nTx,N/nTx]); % 将矩阵转换为[nRx,nTx,N/nTx ]形式
h = 1/sqrt(2)*[randn(nRx,nTx,N/nTx) + 1i*randn(nRx,nTx,N/nTx)]; % 瑞利衰落信道
n = 1/sqrt(2)*[randn(nRx,N/nTx) + 1i*randn(nRx,N/nTx)]; % 0均值高斯白噪声
y = squeeze(sum(h.*sMod,2)) + 10^(-Eb_N0_dB(ii)/20)*n;
End
ML
% 当sHat1 = [1 1] [1 -1] [-1 -1 ] [-1 1],从四组数值中找出最小值及所在位置
sHat1 = [1 1];
sHat1 = repmat(sHat1,[1 ,N/2]);
sHat1Mod = kron(sHat1,ones(nRx,1));
sHat1Mod = reshape(sHat1Mod,[nRx,nTx,N/nTx]); %发送矩阵
zHat1 = squeeze(sum(h.*sHat1Mod,2)) ; %接收矩阵
J11 = sum(abs(y - zHat1),1);%将两行加为一行
sHat2 = [1 -1];
sHat2 = repmat(sHat2,[1 ,N/2]);
sHat2Mod = kron(sHat2,ones(nRx,1));
sHat2Mod = reshape(sHat2Mod,[nRx,nTx,N/nTx]);
zHat2 = squeeze(sum(h.*sHat2Mod,2)) ;
J10 = sum(abs(y - zHat2),1);
sHat3 = [-1 1];
sHat3 = repmat(sHat3,[1 ,N/2]);
sHat3Mod = kron(sHat3,ones(nRx,1));
sHat3Mod = reshape(sHat3Mod,[nRx,nTx,N/nTx]);
zHat3 = squeeze(sum(h.*sHat3Mod,2)) ;
J01 = sum(abs(y - zHat3),1);
sHat4 = [-1 -1];
sHat4 = repmat(sHat4,[1 ,N/2]);
sHat4Mod = kron(sHat4,ones(nRx,1));
sHat4Mod = reshape(sHat4Mod,[nRx,nTx,N/nTx]);
zHat4 = squeeze(sum(h.*sHat4Mod,2)) ;
J00 = sum(abs(y - zHat4),1);
rVec = [J11;J10;J01;J00];
[jj dd] = min(rVec,[],1);
ref = [1 1; 1 0; 0 1; 0 0 ];
ipHat = zeros(1,N);
ipHat(1:2:end) = ref(dd,1);
ipHat(2:2:end) = ref(dd,2);
f = find([ip - ipHat]);%发生错误所在位置
nErrML(ii) = size(find([ip - ipHat]),2);%错误个数ZF算法
% 接收端求伪逆矩阵 G = inv(H^H*H)*H^H ,H^H*H 的矩阵维数为nTx xnTx, 求[a b; c d] 逆 1/(ad-bc)[d -b;-c a]
hCof = zeros(2,2,N/nTx);
hCof(1,1,:) = sum(h(:,2,:).*conj(h(:,2,:)),1); % d
hCof(2,2,:) = sum(h(:,1,:).*conj(h(:,1,:)),1); % a
hCof(2,1,:) = -sum(h(:,2,:).*conj(h(:,1,:)),1); % c
hCof(1,2,:) = -sum(h(:,1,:).*conj(h(:,2,:)),1); % b
hDen = ((hCof(1,1,:).*hCof(2,2,:)) - (hCof(1,2,:).*hCof(2,1,:))); %
ad-bc
hDen = reshape(kron(reshape(hDen,1,N/nTx),ones(2,2)),2,2,N/nTx); %
hInv = hCof./hDen; % inv(H^H*H)
hMod = reshape(conj(h),nRx,N); % H^H
yMod = kron(y,ones(1,2)); % 接收信号矩阵化
yMod = sum(hMod.*yMod,1); % H^H * y
yMod = kron(reshape(yMod,2,N/nTx),ones(1,2)); %
yHat = sum(reshape(hInv,2,N).*yMod,1); % inv(H^H*H)*H^H*y
% 接收端硬判决解码计算统计错误
ipHat = real(yHat)>0;
nErrZF(ii) = size(find([ip- ipHat]),2);
MMSE算法
%W = inv(H^H*H+sigma^2*I)*H^H,H^H*H i维数为[nTx x nTx],求逆1/(ad-bc)[d -b;-c
a]
hCof = zeros(2,2,N/nTx) ;
hCof(1,1,:) = sum(h(:,2,:).*conj(h(:,2,:)),1) + 10^(-Eb_N0_dB(ii)/10); % d
hCof(2,2,:) = sum(h(:,1,:).*conj(h(:,1,:)),1) + 10^(-Eb_N0_dB(ii)/10); % a
hCof(2,1,:) = -sum(h(:,2,:).*conj(h(:,1,:)),1); % c
hCof(1,2,:) = -sum(h(:,1,:).*conj(h(:,2,:)),1); % b
hDen = ((hCof(1,1,:).*hCof(2,2,:)) - (hCof(1,2,:).*hCof(2,1,:))); % ad-bc hDen = reshape(kron(reshape(hDen,1,N/nTx),ones(2,2)),2,2,N/nTx);
hInv = hCof./hDen; % inv(H^H*H)
hMod = reshape(conj(h),nRx,N); % H^H
yMod = kron(y,ones(1,2));
yMod = sum(hMod.*yMod,1); % H^H * y
yMod = kron(reshape(yMod,2,N/nTx),ones(1,2));
yHat = sum(reshape(hInv,2,N).*yMod,1); % inv(H^H*H)*H^H*y
%硬判决及错误
ipHat = real(yHat)>0;
nErrMMSE(ii) = size(find([ip- ipHat]),2);