基于MIMO系统的球形译码算法的matlab性能仿真

目录

一、理论基础

1.1 球形译码概述

1.2 球形译码理论

二、核心程序

三、仿真结论

---

一、理论基础

       在BLAST检测中，目前采用的ZF(迫零) 算法,MMSE(最小均方误差) 算法, OSIC(排序连续干扰抵消) 或ML (最大似然)准则来进行译码。前三种算法,实现起来较简单,但是误码率性能较差;而使用ML检测能得到更好的性能,但是其复杂度较高,不易于实现。基于ML检测的SD(球形译码)算法是一种性能优化,复杂度适中的检测算法。已经证明,采用穷尽搜索的ML 检测算法的复杂度随天线数呈指数增长,而SD算法的复杂度在很大信噪比范围内与天线数呈多项式关系。故SD算法可以用较少的计算量来获得最大似然译码性能。

1.1 球形译码概述

        球面译码（Sphere Decoding）是一种用于多输入多输出（MIMO）通信系统中的信号译码算法，主要用于解码接收到的复杂信号，以恢复发送的信息。球面译码算法旨在在MIMO系统中高效地找到最有可能的发送符号组合，从而降低解码的复杂度。球面译码的核心思想是在复平面上构建一个“球面”，然后通过逐步减小球面半径的方式，将接收到的信号点与可能的发送符号点进行比较，找到最可能的发送符号组合。球面译码不仅可以用于解调QAM、PSK等调制方式，还可以用于解码空时编码（STBC）和空分复用（SFBC）等。

以下是球面译码的主要步骤和概述：

1. 接收信号与发送符号关系建模： 在MIMO系统中，接收到的信号可以表示为接收向量，而发送的符号则是复平面上的点。球面译码算法基于接收信号与可能的发送符号之间的关系进行工作。

2. 球面构建： 在复平面上构建一个球面，初始半径足够大，能够覆盖所有可能的发送符号点。

3. 球面逐步缩小： 球面译码从球面的外围开始逐步缩小，以寻找接收信号最可能对应的发送符号。每次缩小球面半径，都会根据接收信号与可能的发送符号之间的距离，排除掉不太可能的发送符号点。

4. ML（Maximum Likelihood）决策： 在每个缩小的球面上，通过计算接收信号与球面上的各个可能发送符号点之间的距离，找到距离最小的点，即最可能的发送符号。

5. 递归搜索： 在球面逐步缩小的过程中，如果某一步找不到符合条件的发送符号点，就进行递归搜索，扩大搜索半径，以确保能够找到合适的解。

6. 解码： 最后，将找到的最可能的发送符号解码为对应的信息。

        球面译码的优势在于它能够在较低的复杂度下，接近于最大似然解码的性能。然而，球面译码也有一些局限性，例如在高信噪比情况下，可能需要更多的迭代次数来找到合适的解。此外，球面译码在高维空间中（例如大量天线数目的MIMO系统）可能会变得更加复杂。

        尽管球面译码在一些情况下能够取得良好的性能，但在实际应用中，人们也会根据系统复杂度和性能需求考虑其他解码算法，如迭代译码、近似算法等。

1.2 球形译码理论

       球形译码的基本思想是在以一个矢量x 为中心的半径为d的多维球内搜索格点，通过限制或者减少搜索半径从而减少搜索的点数,进而使得计算时间减少。球形译码算法带来的优点在于它不需要象传统的最大似然译码算法那样需要在整个格内对所有的格点进行搜索,而只需要在一个事先设定的有限球形区域进行搜索,如果该区域所包含的点数相对于整个格内的总点数是相当小的,搜索时间就会大大减少。

       影响球形译码的关键问题有：(1) 怎样选择搜索半径d。如果d太大,则球内会包含太多的点,复杂度就会接近或者达到最大似然译码的指数级复杂度。如果d 太小,则球内可能一个格点都不包含,那么球形译码算法将得不到合理的解。(2) 怎样才能判断一个点是否在球内。如果这种判断需要借助每一个格点和矢量之间的距离来判断的话,那么这种方法就不太理想,因为我们需要考察所有的点,所产生的计算量也是指数级的。

       球形译码解决了第2个问题，此处均考虑信号为实数，因为复数可以通过增加一倍的维数，将实部和虚部分开，要判断一个点是否在半径为d的m维球内比较困难。若将m变为1，则从球退化为一个间距，这个点就相当于某根天线发送信号的实部或虚部，这样操作就简单很多，可以知道这个点是否在这个距离内。多根发送天线上的信号的实部和虚部分成很多维，每一维上有可能取值。球形译码算法相当于构建了一棵树，树的第k层节点对应的是落在半径为d，维数为k的球内的格点

 

二、核心程序

function y1 = spheredecodeinftoML(rev,C,H, codebook)
[m,n]=size(H);      %m为Nr，n为Nt
T = zeros(1,n);   

[Q,R] = qr(H);    %QR分解
R = R(1:n, :);    %R= R（1：Nt，：）
ro = Q'*rev;  %y'=Q1'y
for ii = 1:n  %R为上三角矩阵
    if R(ii,ii) == 0
        R(ii,ii) = 0.001;
    end
    for jj = ii+1:n
        R(ii,jj) = R(ii,jj)/R(ii,ii); %为后面的运算做铺垫
    end
    ro(ii) = ro(ii)/R(ii,ii);%相当于ym/Rmm，ym-1/Rm-1...
    R(ii,ii) = R(ii,ii)^2;
end
d = C;    %初始半径
T(n) = 0; %T相当于是dm，m-1，m-2...的值，每维的半径
S = ro;   %相当于ym/Rmm
cn = length(codebook); %码本长度，4qam的话，就为2
len = zeros(1,n);
Y = zeros(n, cn);
ind = zeros(1,n);
y1 = zeros(n,1);
yt = [];
while 1
    k = n;                          %表示从第n维开始计算可能取值
    while k >= 1
       U =sqrt(d/R(k,k))+S(k);%上限值
       ii=cn;
        while ii >= 1 && U < codebook(ii)-10^(-11)
            ii = ii - 1;
        end
       L =-sqrt(d/R(k,k))+S(k);%下限值
       jj=1;
        while jj <= cn && L > codebook(jj)+10^(-11)
            jj = jj + 1;
        end
       len(k) = ii-jj+1;    %第K维的可能值的个数
       Y(k,:) = zeros(1, cn);
       yt = [];
        for kk = 1:len(k)
            yt = [yt codebook(jj+kk-1)];
        end
        [yt1, idex] = sort(abs(yt - S(k)));%按行升序排序
        yt = yt(idex);                     %yt存储信号的可能取值    
        for kk = 1:len(k)
            Y(k,kk) = yt(kk);%从偏差最小的点开始存储，按序，存储信号的可能取值在Y里
        end
        ind(k)=0;%作为标记
        d1 = d+1;
        while d1 > d - 10^(-11)
            ind(k) = ind(k) + 1;
            while ind(k) > len(k)
                if k == n
                    if y1~=zeros(n,1)
                         return       %结束，返回y1    
                    else
                          C=2*C;
                          y1=spheredecodetoML(rev, C, H, codebook);
                          return  
                    end   
                end
                k = k+1;   %表示在第k维中没有找到值，则返回到k+1维中
                ind(k) = ind(k) + 1; %在第k+1维中从下一个点开始搜
            end
            if k > 1
              %  T(k-1) = T(k) - R(k,k)*(S(k)-Y(k,ind(k)))^2;%求dm-1的值，Y(k,ind(k))相当于第k维的取值   采用第二范数
               T(k-1)=max(T(k),R(k,k)*(S(k)-Y(k,ind(k)))^2);                                    % 采用无穷范数
                tmp = 0;
                for ii = k:n
                    tmp = tmp + R(k-1, ii)*Y(ii,ind(ii)) ;
                end
                S(k-1) = ro(k-1) - tmp;              %
                k = k - 1;                           %从第k维降到k-1维
                break;                               %跳出while循环              
            end
          %  d1 = T(n) - T(1) + R(1,1)*(S(1) - Y(1, ind(1)))^2;%k=1时的d1值
         d1 = max(T(1), R(1,1)*(S(1) - Y(1, ind(1)))^2); %k=1时的值    无穷范数
        end %while d1 > d
        if d1 < d-10^(-11)
            for ii=1:n
                y1(ii) = Y(ii, ind(ii));
            end
          d = d1;
        T(n)=0;
        break;
      end
    end %while k>=1
end %while 1 k=n
up212

三、仿真结论