初识MIMO-OFDM（三）：完整的MIMO-OFDM仿真
零.缘起
MIMO-OFDM的仿真是一个大的不能再大的话题了，写得好的仿真平台是可以赚大钱的，我后续也会更新一个著名仿真平台的学习博客，不过我们这里先来学习一下陈老湿在通信家园看到的这个仿真https://zhuanlan.zhihu.com/p/392827532。这个博客主要是希望梳理每一个函数，然后给出2*2Alamouti+OFDM完整的仿真流程。
代码地址：https://github.com/liu-zongxi/MIMO_OFDM_simulation
发射机
一. 发射机步骤一：生成所需要发送的数据
%-----------------------生成比特-----------------%
%-----------------------author:lzx-------------------------%
%-----------------------date:2022年5月5日14点33分-----------------%
function Frame_bit = FrameBitGenerator(N_data, N_user, N_mod, N_symbol)
% 输入
% N_data:一个符号中要发送的符号数
% N_user:用户数
% N_mod:调制数
% N_symbol:一帧中有多少个符号
% 输出
% Frame_bit:一个（1，N_user）的结构体，包含data和num两个field，分别表示每个用户的数据和数据长度
N_bit_per_symbol = N_mod * N_data;
% 可以选择用户比特分布为：
% 1) 所有用户比特相同。
N_bit_per_user_per_symbol = repmat(N_bit_per_symbol/N_user,1,N_user);
% 2) 用户按照比例 (1:N_user)/((1 + N_user)*N_user/2) 发送比特数据
% 目的是为了仿真不同用户的信息比特长度不同的情况,也可以修改得到其他的用户数据比例
% N_bit_per_user_per_symbol = round([1:N_user]/((1 + N_user)*N_user/2) * N_bit_per_symbol);
N_bit_per_user = N_bit_per_user_per_symbol * N_symbol;
for iuser = 1:N_user
% frame_bit_temp{iuser} = rand(N_bit_per_user(iuser), 1) > 0.5;
Frame_bit(iuser).data = rand(N_bit_per_user(iuser), 1) > 0.5;
Frame_bit(iuser).num = N_bit_per_user(iuser);
end
end
该函数主要是生成所需要发送的数据，这里是为了扩展性用了结构体，实际上该仿真根本无法处理不同长度的数据。
一个OFDM符号所承载的比特数是他用于data传输的子载波数*调制数
由于四个用户等分数据，因此每个用户收到数据个数都要被除以N_user
最后我们得到的数据形式是
（N_data*N_mod/N_user*N_symbol, N_user）
二.发射机步骤二：子载波分配
%-----------------------子载波分配-----------------%
%-----------------------author:lzx-------------------------%
%-----------------------date:2022年5月5日16点35分-----------------%
function [index_data_per_user, Frame_zero_padding]=SubcarrierAllocation(Frame_bit_coded, index_data, N_data, N_user,N_symbol, N_mod, type_alloc)
Frame_zero_padding = cell(1,N_user);
index_data_per_user = cell(1,N_user);
% 计算填零的个数
N_subcarrier_per_user = N_data / N_user; % 平均分配,每个用户所被分配的子载波个数
for iuser = 1:N_user
N_zero_padding = N_mod * N_subcarrier_per_user * N_symbol - Frame_bit_coded(iuser).num;
% 子载波分配方式，这里只给出相邻分配
if type_alloc == "neighbour"
index_data_per_user{iuser} = index_data((iuser-1)*N_subcarrier_per_user + 1: iuser*N_subcarrier_per_user)';
end
Frame_zero_padding{iuser} = [Frame_bit_coded(iuser).data; zeros(N_zero_padding, 1)];
end
该函数还是装神弄鬼了，他输出的是针对步骤一所得到的（N_data*N_mod/N_user*N_symbol, N_user）的数据，每一列的每一个symbol使用编号多少的子载波进行发送？
原仿真里给出了多种分配方式，但我这里单纯是实现一下，就是按照顺序分配的
最终输出的index_data_per_user是一个（N_data/N_user，N_user）的矩阵
Frame_zero_padding是一个零填充帧，如果现实情况中不是像我们步骤一中恰好生成那么多比特呢？这时要用0把他填满
三.发射机步骤三：调制
%-----------------------调制函数---------------------%
%-----------------------author:lzx--------------------------%
%-----------------------date:2022年5月5日23点22分-----------------%
function Frame_mod = Modulator(Frame_zero_padding, index_data_per_user, N_data, N_user, N_symbol, N_mod, N_subcarrier)
Frame_mod = zeros(N_subcarrier, N_symbol);
N_subcarrier_per_user = N_data / N_user;
for iuser = 1:N_user
for isymbol = 1:N_symbol
L_symbol = N_subcarrier_per_user* N_mod;
index_symbol = (isymbol-1)*L_symbol+1 : isymbol*L_symbol;
Symbol = Frame_zero_padding{iuser}(index_symbol);
Symbol_premod = reshape(Symbol, N_mod, N_subcarrier_per_user);
Symbol_mod = SymbolModulator(Symbol_premod);
Frame_mod(index_data_per_user{iuser}, isymbol) = Symbol_mod.';
end
end
%-----------------------MIMO-OFDM主函数---------------------%
%-----------------------author:lzx--------------------------%
%-----------------------date:2022年5月5日10点29分-----------------%
function Symbol_mod = SymbolModulator(Symbol_premod)
N_mod = size(Symbol_premod,1);
switch N_mod
% BPSK调制
case 1
matrix_mapping = [-1 1];
Symbol_mod = matrix_mapping(Symbol_premod + 1);
% QPSK调制
case 2
% 比特的映射关系,00:-3/4*pi,01:3/4*pi,10: -1/4*pi,11: 1/4*pi
matrix_mapping = exp(1j*[-3/4*pi 3/4*pi -1/4*pi 1/4*pi]);
index_mod = [2 1]*Symbol_premod;
% 把输入比特映射为符号
Symbol_mod = matrix_mapping(index_mod + 1);
% 8PSK调制
case 3
matrix_mapping = exp(1j*[0 1/4*pi 3/4*pi 1/2*pi -1/4*pi -1/2*pi pi -3/4*pi ]);
index_mod = [4 2 1]*Symbol_premod ;
% 把输入比特映射为符号
Symbol_mod = matrix_mapping(index_mod + 1);
% 16QAM调制
case 4
% 映射关系参见说明文档
m=1;
for k=-3:2:3
for l=-3:2:3
% 对符号能量进行归一化
matrix_mapping(m) = (k+1j*l)/sqrt(10);
m=m+1;
end
end
matrix_mapping = matrix_mapping([0 1 3 2 4 5 7 6 12 13 15 14 8 9 11 10]+1);
index_mod = [8 4 2 1]*Symbol_premod ;
Symbol_mod = matrix_mapping(index_mod + 1);
% 64QAM调制
case 6
% 映射关系参见说明文档
m=1;
for k=-7:2:7
for l=-7:2:7
% 对符号能量进行归一化
matrix_mapping(m) = (k+1j*l)/sqrt(42);
m=m+1;
end
end
matrix_mapping = matrix_mapping(...
[[ 0 1 3 2 7 6 4 5]...
8+[ 0 1 3 2 7 6 4 5]...
24+[ 0 1 3 2 7 6 4 5]...
16+[ 0 1 3 2 7 6 4 5]...
56+[ 0 1 3 2 7 6 4 5]...
48+[ 0 1 3 2 7 6 4 5]...
32+[ 0 1 3 2 7 6 4 5]...
40+[ 0 1 3 2 7 6 4 5]]+1);
index_mod = [32 16 8 4 2 1]*Symbol_premod ;
Symbol_mod = matrix_mapping(index_mod + 1);
end
该函数是一个调制函数
分为两部分，首先在Frame选出一个user一个symbol所需要发送的比特，改变形状后进行调制，调制后放入步骤二所获得的对应的子载波中。
进行调制输入比特为（N_mod, N_subcarrier_per_user），输出一个（1，N_subcarrier_per_user），然后反转后放入到对应子载波中
调制采用的是非常巧妙的方式，首先获得格雷码的映射，然后一一对应，可惜仔细学习一下该函数。
四.步骤四：Alamouti编码
%-----------------------STBC时空块码-----------------%
%-----------------------author:lzx-------------------------%
%-----------------------date:2022年5月5日15点52分-----------------%
function Frame_STBC = STBCCoding(Frame_mod, N_subcarrier, N_symbol, N_Tx)
Frame_STBC = zeros(N_subcarrier, N_symbol, N_Tx);
if (mod(N_symbol,N_Tx))
error('空时编码器输入符号不匹配,子程序st_coding出错');
else
for ispace = 1:N_symbol/N_Tx
if N_Tx == 2
X1=Frame_mod(:,(ispace-1)*N_Tx+1);%取第一列，即第一个OFDM符号
X2=Frame_mod(:,(ispace-1)*N_Tx+2);%取第二列，即第二个OFDM符号
Symbol_STBC = [X1 X2;-conj(X2) conj(X1)];%alamouti编码
elseif N_Tx == 4
X1=Frame_mod(:,(ispace-1)*N_Tx+1);%取第一列，即第一个OFDM符号
X2=Frame_mod(:,(ispace-1)*N_Tx+2);%取第二列，即第二个OFDM符号
X3=Frame_mod(:,(ispace-1)*N_Tx+3);%取第一列，即第一个OFDM符号
X4=Frame_mod(:,(ispace-1)*N_Tx+4);%取第二列，即第二个OFDM符号
Symbol_STBC = [ X1 X2 X3 X4;...
-X2 X1 -X4 X3;...
-X3 X4 X1 -X2;...
-X4 -X3 X2 X1;...
conj(X1) conj(X2) conj(X3) conj(X4);...
-conj(X2) conj(X1) -conj(X4) conj(X3);...
-conj(X3) conj(X4) conj(X1) -conj(X2);...
-conj(X4) -conj(X3) conj(X2) conj(X1)];
end
for iant = 1:N_Tx
Symbol_STBC_per_ant = reshape(Symbol_STBC(:,iant), N_subcarrier, N_Tx);
Frame_STBC(:, (ispace-1)*N_Tx+1:ispace*N_Tx, iant) = Symbol_STBC_per_ant;
end
end
end
这是一个标准的STBC编码，每次取出相邻的两个symbol的符号（同时包含导频和数据，就是一个完整的symbol）（因为这才满足Alamouti编码0时刻和T时刻的要求）。然后进行STBC的编码把原本（N_subcarrier，2）变为（2*N_subcarrier，2）。之后对于编码后的每一列，是由一根天线发送的，把每一列reshape成（N_subcarrier， 2）（这个2是因为2*2编码），放入总矩阵中，当每一列都摆放后，得到一个（N_subcarrier，2，2）的矩阵，这其中第二个2是由于这个2是因为2*2编码，第三个2是由于两根天线，其实是一样的。
当所有的symbol都被取出后，最后得到一个（N_subcarrier，N_symbol，2）的矩阵，这就是编码结果了。
从结果上来看，它将原本（N_subcarrier，N_symbol）个信号生成为了（N_subcarrier，N_symbol，2）个信号，给两个天线分别发送，这就是发射分集。
五.步骤五：OFDM调制
%-----------------------调制函数---------------------%
%-----------------------author:lzx--------------------------%
%-----------------------date:2022年5月6日21点11分-----------------%
function Frame_transmit = OFDMModulator(Frame_pilot, N_sym, N_subcarrier, N_symbol, N_Tx, N_GI)
Frame_transmit = zeros(1,N_sym*N_symbol,N_Tx);
for iant = 1:N_Tx
% ifft乘sqrt(N_subc)以保证变换前后能量不变
% 我们假设频域的样点是在[-fs/2 fs/2]中的, fs是采样频率
% 使用fftshift函数目的是使得变换前的频域样点转换到[0 fs]中,以满足IFFT变换的要求
Frame_TD = sqrt(N_subcarrier) * ifft( fftshift( Frame_pilot(:,:,iant), 1 ) );
Frame_CP = Frame_TD(N_subcarrier - N_GI + 1:N_subcarrier ,:);
Frame_withCP = [Frame_CP; Frame_TD];
% 转换为串行信号
Frame_transmit(:,:,iant) = reshape( Frame_withCP, 1, N_sym*N_symbol);
end
该函数是对信号进行OFDM调制，调制按照天线来进行，将（N_subcarrier，N_symbol）大小的数据在列这一维上做fftshift然后进行ifft，这样还会得到（N_subcarrier，N_symbol）的数据，然后加CP，CP的大小是（N_GI, N_symbol），叠在数据的头上，最后，发送时转换为串行信号，发送一个（1，(N_GI+N_subcarrier)*N_symbol）的信号，多跟天仙的结果是（1，(N_GI+N_subcarrier)*N_symbol，2）.
这里为什么要做并串转换呢？实际上就是按顺序发送symbol，因为天线一次只能发送一个symbol呀~
信道
一.信道步骤一：加噪
这已经是老生常谈了
Power_transmit = var(Frame_transmit);%发送信号功率
N_noise = size(Frame_transmit,2);
noise = NoiseGenerator(EbN0, Power_transmit, N_noise);
Frame_noise = Frame_transmit+noise;
%-----------------------生成高斯噪声---------------------%
%-----------------------author:lzx--------------------------%
%-----------------------date:2022年5月6日22点43分-----------------%
function noise = NoiseGenerator(EbN0, Power_transmit, N_noise)
sigma = sqrt(Power_transmit/(2*EbN0));%标准差sigma
noise_real = randn(1,N_noise,2);
noise_real(:,:,1)=noise_real(:,:,1).*sigma(1);
noise_real(:,:,2)=noise_real(:,:,2).*sigma(2);
noise_imag = randn(1,N_noise,2);
noise_imag(:,:,1)=noise_imag(:,:,1).*sigma(1);
noise_imag(:,:,2)=noise_imag(:,:,2).*sigma(2);
noise=complex(noise_real,noise_imag);%复噪声序列
其实这个代码是比较愚蠢的，生成和发送信号等大的一个噪声矩阵，不详细说了
二.信道步骤二：信道
没有仿真的坐而论道，希望各位大佬指教
这个仿真里。。。。。居然没有信道模型，这其实是我最想知道的一个模型
目前对于信道的仿真，我们是割裂的。
在SISO-OFDM中，我们采用的是一个标准的TDL模型，他的大小是（1，Tau_max），然后进行卷积
在MIMO中，我们使用的是一个平坦信道，他是一个（N_Rx,NTx）的信道，直接进行相乘
那么在MIMO-OFDM中如何进行结合？因为每一个子载波都应该对应着一个子信道，信道的数量应该是(N_subcarrier, N_Tx*N_Rx);但这应该如何仿真呢？
接收机
首先要明确，接收机有N_user个
一. 步骤一：OFDM解调
%-----------------------OFDM的接收---------------------%
%-----------------------author:lzx-------------------------%
%-----------------------date:2022年5月5日16点35分-----------------%
function Frame_recieve = OFDMDemodulator(Frame_noise, N_sym, N_subcarrier,N_symbol,N_Rx, N_GI)
Frame_recieve = zeros(N_subcarrier,N_symbol,N_Rx);
for iant = 1:N_Rx
Frame_symbol = reshape(Frame_noise(1,:,iant), N_sym, N_symbol);
Frame_noGI = Frame_symbol(N_GI+1:end, :);
% fft乘1/sqrt(N_subc)以保证变换前后能量不变
% 我们假设频域的样点是在[-fs/2 fs/2]中的, fs是采样频率
% fftshift目的是使得变换后的频域样点在[-fs/2 fs/2]中,而不是[0 fs]中
Frame_recieve(:,:,iant) = fftshift(1/sqrt(N_subcarrier) * fft( Frame_noGI ), 1);
end
% Frame_recieve = Frame_FD(:, 1:N_symbol ,:); % 数据OFDM符号, 包括导频符号
作为天线，自然应该按天线处理
接收到的信号重新串并转换为（N_sym, N_symbol)，因为我们肯定是一个symbol一个symbol接收到的，去掉CP然后反做fft即可
最后得到的输出是（N_subcarrier，N_symbol），按天线接收后，得到的是（N_subcarrier，N_symbol，2）不过我这里存在疑问，接收端如何能分开发射端两根天线？我对此感到不解
二. 接收机二：STBC解码
%-----------------------STBC解码---------------------%
%-----------------------author:lzx-------------------------%
%-----------------------date:2022年5月7日10点07分-----------------%
function Frame_decoded = STBCDecoding(Frame_recieve, H, N_subcarrier, N_Tx,N_Rx, N_symbol)
Frame_decoded = zeros(N_subcarrier, N_symbol);
% 把H_freq转化为空时译码器的输入格式,为一个N_subc*N_ant_pair的矩阵,每列表示:
% 1-->1 ,1-->2,...,1-->N_Rx_ant, ... ,N_Tx_ant-->1, N_Tx_ant-->2,..., N_Tx_ant-->N_Rx_ant
if (N_Tx == 2)&&(N_Rx == 2)
for ispace = 1:N_symbol/N_Tx
% 构造输入进空时译码器的符号,用t表示时间号,a表示天线号,其格式为:
% [Recv(t1,a1) Recv(t2,a1) Recv(t1,a2) Recv(t2,a2)].
R = [];
for iant = 1:N_Rx
R = [R Frame_recieve(:,(ispace-1)*N_Tx+1:ispace*N_Tx,iant) ];%对于2X2MIMO来说，R矩阵即为[Xe -conj(Xo) Xo conj(Xe)]
end
H11=H(:,1);%1发送天线--》1接收天线信道参数
H12=H(:,2);%1发送天线--》2接收天线信道参数
H21=H(:,3);%2发送天线--》1接收天线信道参数
H22=H(:,4);%2发送天线--》2接收天线信道参数
R11=R(:,1);%1发送天线发送的第一个符号
R12=R(:,2);%1发送天线发送的第二个符号
R21=R(:,3);%2发送天线发送的第一个符号
R22=R(:,4);%2发送天线发送的第二个符号
for i=1:1:N_subcarrier
X1(i,1)=(R11(i)*conj(H11(i))+conj(R12(i))*H21(i)+R21(i)*conj(H12(i))+conj(R22(i))*H22(i))/( H11(i)*conj(H11(i)) + H21(i)*conj(H21(i)) + H12(i)*conj(H12(i)) + H22(i)*conj(H22(i)));
X2(i,1)=(R11(i)*conj(H21(i))-conj(R12(i))*H11(i)+R21(i)*conj(H22(i))-conj(R22(i))*H12(i))/( H11(i)*conj(H11(i)) + H21(i)*conj(H21(i)) + H12(i)*conj(H12(i)) + H22(i)*conj(H22(i)));
end
Sybmol_decoded = 2*[X1 X2]; % 非常扯淡，因为上面的分母是错误的
Frame_decoded(:,(ispace-1)*N_Tx+1:ispace*N_Tx) = Sybmol_decoded;
end
end
首先我要申明，我认为该仿真是错的Sybmol_decoded = 2*[X1 X2]; % 非常扯淡，因为上面的分母是错误的为什么莫名其妙的×2呢，因为上面( H11(i)*conj(H11(i)) + H21(i)*conj(H21(i)) + H12(i)*conj(H12(i)) + H22(i)*conj(H22(i)));是错误的，应该只选取其中两个。代码的实现更是丑陋，所以我也没有重写的欲望了。她也没有解决我想知道的信道的问题，
对于每两列进行一个操作，再加上两根天线，生成会原本的（N_subcarrier*2,2）的矩阵，然后解算。的总之他会把原本（N_subcarrier，N_symbol，2）的信号接收后得到N_subcarrier，N_symbol).
三.接收机步骤三：解调
%-----------------------QAM解调-----------------%
%-----------------------author:lzx-------------------------%
%-----------------------date:2022年5月7日11点51分-----------------%
function Frame_demod = Demodulator(Frame_decoded, index_data_iuser, N_mod, N_symbol)
Frame_demod = [];
for isymbol = 1:N_symbol
Symbol_demod = SymbolDemodulator(Frame_decoded(index_data_iuser, isymbol).', N_mod);
Frame_demod = [Frame_demod; Symbol_demod(:)];
end
%-----------------------每个符号的QAM解调-----------------%
%-----------------------author:lzx-------------------------%
%-----------------------date:2022年5月7日11点51分-----------------%
function Symbol_demod = SymbolDemodulator(Symbol, N_mod)
switch N_mod
% BPSK解调
case 1
Symbol_demod = real(Symbol) > 0;
% QPSK解调
case 2
% 由QPSK的星座图可以观察到
bit0 = real(Symbol) ;
bit1 = imag(Symbol) ;
% 得到2行, 列数为符号数的输出矩阵
Symbol_demod(1,:) = bit0 > 0;
Symbol_demod(2,:) = bit1 > 0;
% 8PSK解调
case 3
% 参见8PSK的星座图
bit0 = -imag( Symbol * exp(1j*pi/8)) ;
% bit1和bit2解调,都需要进行星座旋转
bit1 = -real(Symbol * exp(1j*pi/8)) ;
bit2 = [];
for isymbol = 1:length(Symbol)
tmp = Symbol(isymbol) * exp(-1j*pi/8);
if ((real(tmp) <0) && (imag(tmp) >0)) || ((real(tmp) >0) && (imag(tmp) <0))
bit2 = [bit2 0];
else
bit2 = [bit2 1];
end
end
Symbol_demod(1,:) = bit0 >0;
Symbol_demod(2,:) = bit1 >0;
Symbol_demod(3,:) = bit2 ; % 已经硬判决
% 16QAM解调
case 4
bit0 = real(Symbol);
bit2 = imag(Symbol);
% 以bit1的生成来说明方法:
% 2/sqrt(10) 为临界值, abs(real(sym))大于此, 则bit1为负,硬判决得到0 ; 反之为正
bit1 = 2/sqrt(10)-(abs(real(Symbol)));
bit3 = 2/sqrt(10)-(abs(imag(Symbol)));
Symbol_demod(1,:) = bit0 > 0;
Symbol_demod(2,:) = bit1 > 0;
Symbol_demod(3,:) = bit2 > 0;
Symbol_demod(4,:) = bit3 > 0;
% 64QAM解调
case 6
bit0 = real(Symbol);
bit3 = imag(Symbol);
bit1 = 4/sqrt(42)-abs(real(Symbol));
bit4 = 4/sqrt(42)-abs(imag(Symbol));
for m=1:size(Symbol,2)
for k=1:size(Symbol,1)
if abs(4/sqrt(42)-abs(real(Symbol(k,m)))) <= 2/sqrt(42)
bit2(k,m) = 2/sqrt(42) - abs(4/sqrt(42)-abs(real(Symbol(k,m))));
elseif abs(real(Symbol(k,m))) <= 2/sqrt(42)
bit2(k,m) = -2/sqrt(42) + abs(real(Symbol(k,m)));
else
bit2(k,m) = 6/sqrt(42)-abs(real(Symbol(k,m)));
end
if abs(4/sqrt(42)-abs(imag(Symbol(k,m)))) <= 2/sqrt(42)
bit5(k,m) = 2/sqrt(42) - abs(4/sqrt(42)-abs(imag(Symbol(k,m))));
elseif abs(imag(Symbol(k,m))) <= 2/sqrt(42)
bit5(k,m) = -2/sqrt(42) + abs(imag(Symbol(k,m)));
else
bit5(k,m) = 6/sqrt(42)-abs(imag(Symbol(k,m)));
end
end
end
Symbol_demod(1,:) = bit0 > 0;
Symbol_demod(2,:) = bit1 > 0;
Symbol_demod(3,:) = bit2 > 0;
Symbol_demod(4,:) = bit3 > 0;
Symbol_demod(5,:) = bit4 > 0;
Symbol_demod(6,:) = bit5 > 0;
end
这里对于QAM，我没有自己写过就不班门弄斧了，但我大致知道他的思路
解调时，也是对于每一个symbol的每一个用户的那一部分子载波进行解调，挑选出属于该接收机 子载波然后全部并联起来，最后得到一个（N_data*Nysmbol, 1）便于计算误码率。文章来源：https://www.toymoban.com/news/detail-402130.html
总结
总的来说，这个仿真我觉得不行!而且我认为OFDM的发送解调和STBC的发送解调顺序可能有问题。如果后续看到更好的仿真，我会继续更新的！文章来源地址https://www.toymoban.com/news/detail-402130.html
到了这里，关于初识MIMO-OFDM（三）：完整的MIMO-OFDM仿真的文章就介绍完了。如果您还想了解更多内容，请在右上角搜索TOY模板网以前的文章或继续浏览下面的相关文章，希望大家以后多多支持TOY模板网！
