第7章：OFDM 信道估计与均衡（2）

本文对应的完整可运行代码下载地址：123kevin456/OFDM-。

一、瑞利衰落信道介绍

第（1）内容《第7章：OFDM 信道估计与均衡（1）》讲到OFDM经过高斯白噪声信道下的误码率分析，这一节重点关注OFDM经过多径衰落信道的误码率分析。从这里也将看到保护间隔、信道均衡对于接收端正确解调的决定性作用。

我前面已经多次提到均衡，那么均衡（equalization）是什么呢？均衡就是把符号间干扰去除的过程。

（那么符号间为什么会有干扰呢？这个我在《第1章：BPSK调制解调器仿真》解释过，不清楚的同学可以去前面复习一下。）

继续介绍瑞利衰落信道，当然会涉及本科的《无线通信》课程了。

无线信道中，发射信号会可能有直射、折射、反射路径到达接收端，由于不同路径的长度不一样，所以导致同一信号的不同分量到达接收端的时间不一样。

那么问题来了，为什么会有直射、反射等路径呢？这是因为天线指向性，有主瓣、旁瓣之分。

举个通俗例子来理解，如果将主瓣指向接收机，旁瓣的能量便不会指着接收机，朝别的方向发射去了。

关于天线的概念，发送天线增益、接收天线增益、方向图的概念有没有？本科的《微波技术与天线》都会讲到。

明白了为什么有多径，当然，肯定要有发射体、折射体，不然旁瓣打出去的能量就白白浪费了，没有办法到达接收端了。

多径会对接收端造成什么影响呢？这是一个关键性问题。

图1 多径衰落信道

在仿真过程中，不妨以第一个路径，即最短路径，到达接收端的时间和功率衰减作为基准，第二、第三等路径相比第一路径做相应变化。

举个数字例子帮助理解：

PowerdB=[0 -8 -17 -21 -25]; % 信道抽头功率特性
Delay=[0 3 5 6 8];          % 信道时延,示例

由上图的PowerdB和Delay，即第二径相比第一径晚到3s，第三径相比第一径晚到5s，以此类推。功率上来说，第二径相比第一径衰减8dB，第三径相比第一径衰减17dB，以此类推。

你可能会问了，为什么要采取这种相对第一径的方式呢？好用呗，哈哈哈哈。

PowerdB和Delay这两个参数，将直接应用到后面的仿真代码中。

先看到衰落信道的特性，以及如何产生：

图2 衰落信道的统计特性

上面三张截图来自《MIMO-OFDM无线通信技术及MATLAB实现》（这本书会经常被我点名），每条路径均考虑了信道增益、时延和多普勒频移。在下面的代码实现中，暂不考虑多普勒频移。

莱斯衰落和瑞丽衰落的概念要有，这是两种经典的多径衰落。

在时域上，观测信号是输入信号和系统响应的线卷积，表达为： $y=Hx+n$ 。

其中H是线卷积矩阵，长下面这样：

接收端面对这样的信道条件，有两种处理方式：

（1）从时域角度看

接收端若可以估计不同路径的信道增益，便可以将接收到的信号乘以复数共轭的冲击响应。当然，这要求接收端能够存下大量的连续数据用来做均衡。

（2）从频域角度看

发射信号经过不同路径到达接收端，有可能会起到增强作用、也有可能会起到衰减作用。对于窄带信号来说，当信号带宽大于信道的相干带宽时，信号便会发生严重的失真(distortion)。而当信号带宽小于信道的相干带宽时，信号发生失真便很小。

考虑到宽带移动的场景，常常需要有较高的数据传输速率，结合上面讲的多径衰落中信号带宽与信道带宽的关系，于是将高速的数据传输被转化成多通道的并行传输。

不同的通道数据如何来做区分呢？

数据之间采用多路复用技术来区分不同子信道，比如可以用频率、用码来区分。（关于多路复用和多址接入的区别，这两个概念不是一个意思，思考几秒钟？然后自行百度？）

上面这几段话是什么意思呢？接下来我将举例说明。

对于给定速率要求的系统来说，速率越高，意味着每个符号占据的时间越短，占用的带宽便会很多。假设你想以1MBaud的速率发送符号，而现在只有一个通道供你发射，那么单个符号持续时间是1us，占据带宽假设是0.5MHz。

倘若有2个通道供你发射，每个通道便可以按0.5MBaud速率发送，那么单个符号持续时间此时是2us，只需要各占0.25MHz。

那你可能会说，2个信道也是需要0.5MHz发送呀。

结合下面这个图来看，确实是这样。但是并行传输能将频率选择性衰落，转化成多个平坦衰落，这便是OFDM的优势了。

图3 多径衰落信道冲击响应的时域和频域图

二、OFDM经过多径衰落信道的误码率分析

接下来直接上代码，结合代码来讲解OFDM经过多径信道后的误码率情况，以及均衡是如何发挥作用的。

%%%%%%%%%%%%%%%%%%%%% OFDM仿真 %%%%%%%%%
%%%%%%%%%%%%%%%%%%%%% ofdm_fading_sim2.m  %%%%%%%%%
%%%%%%%%%  data:2020年10月16日  author:飞蓬大将军 %%%%%%%%%%


%%%%%程序说明
%%%多径衰落信道下的OFDM传输
%%%调制方式：QPSK
%%%编码方式：无
%%%接收端不做捕获和同步，也不做信道估计，假设知道信道条件，后续版本考虑信道估计

%%%%    仿真环境
%软件版本：MATLAB R2019a

%********************** 程序主体 ************%

%%%%%%%%%%%%%%%%%%%%%  参数设置   %%%%%%%%%%%%%%%%%%%
para = 128;   %Number of parallel channel to transmit
fftlen = 128;  %FFT length
noc = 128;    %Number of carrier 
nd = 6;  %Number of information OFDM symbol for one loop
ml = 2;   %Modulation：QPSK
sr = 250000;  %Symbol rate 符号速率
br = sr.*ml;  %Bit rate per carrier
gilen = 32; %length of guard interval 
% gilen = 0; %length of guard interval 
ebn0_temp = 0:2:20;
ber_fading = zeros(1,length(ebn0_temp));
for kkk = 1:length(ebn0_temp)
    
ebn0 = ebn0_temp(kkk);  %Eb/No 


%%%%%%%%%%%%%%%%%%%%%%%%%%  Fading initialization %%%%%%%%%%%%%%

PowerdB=[0 -8 -17 -21 -25]; % 信道抽头功率特性
Delay=[0 3 5 6 8];          % 信道时延,示例
% Delay=[0 3 5 56 78];          % 信道时延
Power=10.^(PowerdB/10);     % 信道抽头功率特性 '线性'
Ntap=length(PowerdB);       % 信道抽头数
Lch=Delay(end)+1;           % 信道长度


%%%%%%%%%%%%%%%%%%%%%  主循环 %%%%%%%%%%%%%

nloop = 10000; %Number of sumulation loops 
noe = 0;   %Number of error data  
nod = 0;   %Number of transmitted data 
eop = 0;   %Number of error packet  
nop = 0;   %Number of transmitted packet


for iii = 1:nloop
    %%%%%%%%%%%%%%%%%  发射机  %%%%%%%%%%%%%%%%%%%
   seldata = rand(1,para*nd*ml)>0.5;  %串行数据

%     seldata = ones(1,para*nd*ml); 
    paradata = reshape(seldata,para,nd*ml); %串并转换

    [ich,qch] = qpskmod(paradata,para,nd,ml); %调制
    kmod = 1/sqrt(2);
    ich1 = ich.*kmod;
    qch1 = qch.*kmod;
    %%%%%%%%%%%% IFFT %%%%%%%%%%
    x = ich1 + qch1 *1j;
    
%     spow1 = sum(sum(ich1.^2+ qch1.^2))/nd./para;
    
    y = ifft(x);
    ich2 = real(y);
    qch2 = imag(y);
    
%     spow2 = sum(sum(ich2.^2+ qch2.^2))/nd./para;
    
    
    %%%%%%%% 添加保护间隔 %%%%%%%%
    [ich3,qch3] = giins(ich2,qch2,fftlen,gilen,nd);
    fftlen2 = fftlen + gilen;
    
    
    %****************  Attenuation Calculation **************
    %%%方式一：
    spow = sum(ich3.^2+ qch3.^2)/nd./para;
    attn = 0.5*spow*sr/br*10.^(-ebn0/10);
    attn = sqrt(attn);
    
    %%%方式二：
%     snr = ebn0 + 10*log10(2);
%     attn = sqrt(10.^(-snr/10)*spow/2);
    %以上两种方式表达是一样的
    
    
    %***************  衰落信道 Fading channel ***************%
    channel = (randn(1,Ntap) + 1j * randn(1,Ntap)).*sqrt(Power/2);
    h = zeros(1,Lch);
    h(Delay+1) = channel;
    y = conv(ich3 + 1j*qch3,h);
    ifade = real(y(:,1:length(ich3)));
    qfade = imag(y(:,1:length(ich3)));
    
    %***********************  接收机 *******************%
    %%%%%%%%%% AWGN addition %%%%%%%%%
    [ich4,qch4] = comb(ifade,qfade,attn);

    %%%%%%%% 去掉保护间隔 %%%%%%%%
    [ich5,qch5] = girem(ich4,qch4,fftlen2,gilen,nd);
    
    
    %%%%%%%%%%%%%%  FFT %%%%%%%%
    rx = ich5 + qch5.*1i;
    ry = fft(rx);
%     ich6 = real(ry);
%     qch6 = imag(ry);
    
    %%%%%%%%  信道均衡 %%%%
    %%%注意A的共轭转置和转置的区别，前者是A'，后者是A.'
    H = fft([h,zeros(1,fftlen-Lch)].');
    
    for number = 1:nd
        ch6(:,number) = ry(:,number)./H; 
    end
    
    ich6 = real(ch6);
    qch6 = imag(ch6);
    
    %%%%%%%%%%%%%% demoluation %%%%%%%%%%%%%%
    ich7 = ich6./kmod;
    qch7 = qch6./kmod;
    demodata = qpskdemod(ich7,qch7,para,nd,ml);
    
    %%%%%%%%%%%%%% 并串转换 %%%%%%%%%
    demodata1 = reshape(demodata,1,para*nd*ml);
    
    %%%%%%%%%%%%%%% Bit Error Rate %%%%%%%%%%%
    noe2 = sum(abs(demodata1-seldata));
    nod2 = length(seldata);
    
    %%%%cumulative the number of error and data in noe and nod
    noe = noe + noe2;
    nod = nod + nod2;
     
     %%%计算PER
    if noe2~=0
        eop = eop +1;

        
    end

    nop = nop + 1;
%     fprintf('%f\t%e\t%d\n',iii,noe2/nod2,eop);
    
%      fprintf('%f\t%e\t%d\n',iii,noe2/nod2,eop);
%     
    
end
%**************   output result *************
per = eop/nop;
ber = noe/nod;

ber_fading(1,kkk) = ber;
% fprintf('%f\t%e\t%e\t%d\t\n',ebn0,ber,per,nloop);
% fid = fopen('BERofdm.dat','a');
% fprintf(fid,'%f\t%e\t%e\t%d\t\n',ebn0,ber,per,nloop);
% fclose(fid);


end

%************************* 画误码率曲线进行对比 ******************%
rayleign_one_path_theory = ber_temp(ebn0_temp+1); 
semilogy(ebn0_temp,rayleign_one_path_theory,'-*',ebn0_temp,ber_fading,'-+');
xlabel('比特信噪比');
ylabel('误码率');
title('多径衰落信道下误码率仿真曲线');
legend('理论曲线','实验曲线');
grid on;
%***********   多径瑞利衰落下的OFDM采用BPSK调制的误码率值 **********%

%%%%%%%%%%%%%%%%%%%%%%%      理论值          **************%
%%%%%%%%%%%%%     EbN0(dB)      误码率        
%%%%%%%%%%%%%       3        0.125000000000000
%%%%%%%%%%%%%       4        0.100000000000000
%%%%%%%%%%%%%       5        0.0833333333333333
%%%%%%%%%%%%%       6        0.0714285714285715
%%%%%%%%%%%%%       7        0.0625000000000000
%%%%%%%%%%%%%       8        0.0555555555555556
%%%%%%%%%%%%%       9        0.0500000000000000
%%%%%%%%%%%%%      10        0.0454545454545455

%%%%%%%%%%%%%%%%%   结论    %%%%%%%%%%%%%%%%%%%%
%完成了OFDM经过多径衰落信道的误码率仿真
%OFDM中引入保护间隔，是一种冗余信息，因此相比于理论误码率曲线有10*log（160/128）=0.969dB的损失
%信道均衡中，已经将信道带来的影响补偿上
%2020年11月11日

运行上面的代码，可以得到以下结果：

图4 多径衰落信道理论误码率和实验误码率对比图

关于上面的代码有以下几点需要注意：

（1）瑞利衰落信道的理论误码率表达式

$P_e=\frac{M-1}{M{\log }_{2}M}\left(1-\sqrt{\frac{3\gamma {\log }_{2}M/\left(M^2-1\right)}{3\gamma {\log }_{2}M/\left(M^2-1\right)+1}}\right)$

其中，M为调制阶数， $\gamma =\frac{E_b}{N_0}$ ，在类似这样的误码率表达式中， $\frac{E_b}{N_0}$ 都是十进制的数，因此在仿真时常常有将dB值转化为线性值这一步。

（2）上述信道均衡的理论基础是在上一篇讲过的循环卷积和线性卷积概念，均衡也有时域均衡和频域均衡之分，这次采用的是频域均衡。

（3）假设接收端已经估计出来各条路径的信道增益，后续再考虑插入导频以及怎么估计出来的。

有了上面的代码，还可以做以下工作：

（1）如果接收端不做均衡的话，可以发现在不同信噪比条件下的误码率一直都是最大，均在0.5附近。（这个我试过了，你也可以自己试一下，嘿嘿）

（2）可以改变各个路径到达接收端的时间与CP长度的相对大小关系，来观察误码率情况。比如CP=0，可得下图：

图5 CP长度为0对误码率的影响

还可以调整PowerdB和Delay的值，来观察误码率情况。

PowerdB=[0 -8 -17 -21 -25]; % 信道抽头功率特性
Delay=[0 3 5 6 8];          % 信道时延,示例

三、总结

前面我推荐过《Simulation and Software Radio for Mobile Communications》，里面让我最受触动的是作者的序言，内容如下：

这也再一次体现，经典书籍为何被称为“经典”了。

上面的代码中，还没有考虑**信道估计、以及多普勒频移。关于OFDM插入导频来做信道估计，可以时域插入、也可以频域插入，我下次再讲。

也欢迎读者朋友就相关技术问题与我交流，一起学习，共同进步。请你也别忘了把这篇文章分享给你身边正在学习通信专业的同学们，也许能够帮到Ta。这是《陈老湿·通信MATLAB》仿真的第7章，期待下次更新见！