ICode9
精准搜索请尝试： 精确搜索
首页
编程语言>
数据库
系统相关
互联网
其他分享
Python
Java
JavaScript
android
PHP
首页 > 其他分享> 文章详细

VLS系统-16QAM-OFDM系统仿真实现（MATLAB）

2021-09-25 17:02:06  阅读：285  来源： 互联网

标签：received OFDM carr TS length MATLAB 16QAM CP Nt


参考：基于16QAM调制的OFDM系统仿真实现（MATLAB）_hyl1181的博客-CSDN博客

参考的文章代码缺少QAM调制解调函数，不可以运行，经调试可成功运行。代码请见下文

系统图如下：

 

 

 

%%  https://blog.csdn.net/hyl1181/article/details/110942407
%% 已调试 可成功运行
clc
clear all;
close all;
Nt_carr=256;      %子载波数=FFT点数---256
Np_carr=Nt_carr/2-1; %实际子载波数---127
Sig_per_carr=500; %每子载波含符号数---500
M = 16;
K=log2(M);      %每符号含比特数,16QAM调制---4
CP_Ratio=1/8;
CP_length=CP_Ratio*Nt_carr;%循环前缀长度
TS_length=Nt_carr+CP_length;%训练序列长度
SNR=20;

%==================================================
%================信号产生===================================
baseband_out_length=Np_carr * Sig_per_carr * K;  %所输入的比特数目  500*256*4
% rand( 'twister',0);
baseband_out=round(rand(1,baseband_out_length));%输出待调制的二进制比特流

%==============16QAM调制====================================

      x_reshape=reshape(baseband_out,K,length(baseband_out)/K)';                          %对数据流进行分组，对于16QAM，则每4位一组
      x_de=bi2de(x_reshape,'left-msb');                             %转化成10进制，作为qammod的输入
      x_symbol=reshape(x_de,Np_carr,Sig_per_carr);
     %% ---------M-QAM调制--------
       complex_carrier_matrix=qammod(x_de,M,'Gray','UnitAveragePower',true);  
       d=reshape(complex_carrier_matrix,Np_carr,Sig_per_carr); 
figure;
plot(complex_carrier_matrix,'*r');%16QAM调制后星座图
title('16QAM调制星座图');
axis([-5,5,-5,5]);
grid on   %显示网格线
axis square

%==============串并变换====================================
%矩阵转置时附加有取共轭
complex_carrier_matrix1=reshape(complex_carrier_matrix',Np_carr,Sig_per_carr)';%串并变换Sig_per_carr*Nt_carr 矩阵
%complex_carrier_matrix1=conj(reshape(complex_carrier_matrix',Nt_carr,Sig_per_carr)')%两次转置

%==============埃尔米特映射====================================
carriers=(1:Np_carr)+1;%共轭对称子载波映射  复数数据对应的IFFT点坐标
conjugate_carriers=Nt_carr-carriers+2;%共轭对称子载波映射  共轭复数对应的IFFT点坐标
IFFT_modulation=zeros(Sig_per_carr,Nt_carr);%添0组成IFFT_bin_length IFFT 运算
IFFT_modulation(:,carriers )=complex_carrier_matrix1 ;%未添加导频信号 ，子载波映射在此处
IFFT_modulation(:,conjugate_carriers )=conj(complex_carrier_matrix1);%共轭复数映射

%=================IFFT===========================
time_wave_matrix=ifft(IFFT_modulation,Nt_carr,2);%OFDM调制 即IFFT行变换
%时域波形矩阵，行为每载波所含符号数，列IFFT点数，N个子载波映射在其内，每一行即为一个OFDM符号

%=====================计算峰均比PARR====================================
PAPR=10*log10(Sig_per_carr*Nt_carr*(max(max(time_wave_matrix)))^2/sum(sum(time_wave_matrix.^2)));

%==============降PAPR预操作====================================
b=[exp(pi/4*1i) exp(5*pi/4*1i)];
for i=3:Np_carr
    complex_carrier_matrix1(:,i)=b(1)*complex_carrier_matrix1(:,i);
end
PAPR_temp=zeros(1,4);
for i=1:2
    for j=1:2
        complex_carrier_matrix1_temp=[b(i)*complex_carrier_matrix1(:,1) b(j)*complex_carrier_matrix1(:,2) complex_carrier_matrix1(:,3:Np_carr)];
        carriers=(1:Np_carr)+1;
        conjugate_carriers=Nt_carr-carriers+2;
        IFFT_modulation=zeros(Sig_per_carr,Nt_carr);
        IFFT_modulation(:,carriers )=complex_carrier_matrix1_temp ;
        IFFT_modulation(:,conjugate_carriers )=conj(complex_carrier_matrix1_temp);
        time_wave_matrix=ifft(IFFT_modulation,Nt_carr,2);
        PAPR_temp((i-1)*2+j)=10*log10(Sig_per_carr*Nt_carr*(max(max(time_wave_matrix)))^2/sum(sum(time_wave_matrix.^2)));
    end
end
[~,index]=min(PAPR_temp);
complex_carrier_matrix1_temp=[b(ceil(index/2))*complex_carrier_matrix1(:,1) b(index-2*(ceil(index/2)-1))*complex_carrier_matrix1(:,2) complex_carrier_matrix1(:,3:Np_carr)];
carriers=(1:Np_carr)+1;
conjugate_carriers=Nt_carr-carriers+2;
IFFT_modulation=zeros(Sig_per_carr,Nt_carr);
IFFT_modulation(:,carriers )=complex_carrier_matrix1_temp ;
IFFT_modulation(:,conjugate_carriers )=conj(complex_carrier_matrix1_temp);
time_wave_matrix=ifft(IFFT_modulation,Nt_carr,2);
        
figure;
stem(0:Nt_carr-1, abs(IFFT_modulation(2,1:Nt_carr)),'b*-')%第一个OFDM符号的频谱
grid on
axis ([0 Nt_carr -0.5 4.5]);
ylabel('Magnitude');
xlabel('IFFT Bin');
title('第一个OFDM符号各频率对应的幅度（幅频响应）');

figure;
stem(0:Nt_carr-1, (180/pi)*angle(IFFT_modulation(2,1:Nt_carr)), 'b*-')
hold on
plot(0:Nt_carr-1, (180/pi)*angle(IFFT_modulation(2,1:Nt_carr)), 'go')
axis ([0 Nt_carr -200 +200])
grid on
ylabel('Phase (degrees)')
xlabel('IFFT Bin')
title('第一个OFDM符号各频率对应的相位（相频响应）')

figure;
subplot(211)
plot(0:Nt_carr-1,time_wave_matrix(1,:));%第一个OFDM符号的时域波形
axis([0,Nt_carr-1+CP_length, -0.6, 0.6]);
grid on;
ylabel('Amplitude');
xlabel('Time');
title('发送端第一个OFDM符号的时域波形');

%=====================添加循环前缀CP====================================
CP=time_wave_matrix(:,Nt_carr-CP_length+(1:CP_length));
time_wave_matrix_add_CP=[CP,time_wave_matrix];
subplot(212)
plot(0:Nt_carr-1+CP_length,time_wave_matrix_add_CP(1,:));%第一个OFDM符号的时域波形
axis([0,Nt_carr-1+CP_length, -0.6, 0.6]);
grid on;
ylabel('Amplitude');
xlabel('Time');
title('发送端第一个OFDM符号的时域波形（加CP）');

%=========================并串转换======================================
time_wave_sequence=reshape(time_wave_matrix_add_CP',(Nt_carr+CP_length)*Sig_per_carr,1)';

%=====================添加训练序列TS====================================
% TS_0=zeros(1,TS_length);
% BPSKTable=[-1-8i,1+8i];
% TS_0=2*sqrt(65)*randi([0,1],1,TS_length);
% TS_0=5*sqrt(2)*randi([0,1],1,TS_length);%干扰序列平均功率为TS序列的五倍
% TS_0=2*randi([0,1],1,TS_length)-1;%干扰序列平均功率与TS序列相当
TS_0=sqrt(2)*randi([0,1],1,TS_length);
BPSKTable=[-1,1];
TS_Schmidl_1=BPSKTable(randi([0,1],1,Nt_carr/2)+1);
TS_Schmidl=[TS_Schmidl_1(Nt_carr/2-CP_length+1:Nt_carr/2) TS_Schmidl_1 TS_Schmidl_1];
TS_Minn_1=BPSKTable(randi([0,1],1,Nt_carr/4)+1);
TS_Minn=[-TS_Minn_1(Nt_carr/4-CP_length+1:Nt_carr/4) TS_Minn_1 TS_Minn_1 -TS_Minn_1 -TS_Minn_1];
TS_Park_1=BPSKTable(randi([0,1],1,Nt_carr/4)+1);
TS_Park_2=zeros(1,Nt_carr/4);
for i=1:Nt_carr/4
    TS_Park_2(i)=conj(TS_Park_1(Nt_carr/4-i+1));
end
TS_Park=[conj(TS_Park_2(Nt_carr/4-CP_length+1:Nt_carr/4)) TS_Park_1 TS_Park_2 conj(TS_Park_1) conj(TS_Park_2)];
TS=[TS_0 TS_Schmidl TS_0 TS_Minn TS_0 TS_Park TS_0];
time_wave_sequence_add_TS=[TS time_wave_sequence];

figure;
subplot(311)
plot(0:Nt_carr-1+CP_length,TS_Schmidl);%训练序列波形
axis([0,Nt_carr-1+CP_length, -1.5, 1.5]);
grid on;
ylabel('Amplitude');
xlabel('Time');
title('Schmidl训练序列波形');

subplot(312)
plot(0:Nt_carr-1+CP_length,TS_Minn);%训练序列波形
axis([0,Nt_carr-1+CP_length, -1.5, 1.5]);
grid on;
ylabel('Amplitude');
xlabel('Time');
title('Minn训练序列波形');

subplot(313)
plot(0:Nt_carr-1+CP_length,TS_Park);%训练序列波形
axis([0,Nt_carr-1+CP_length, -1.5, 1.5]);
grid on;
ylabel('Amplitude');
xlabel('Time');
title('Park训练序列波形');

%=====================经过高斯噪声信道====================================
figure;
received_time_wave_sequence_add_TS=awgn(time_wave_sequence_add_TS,SNR,'measured');
received_TS=received_time_wave_sequence_add_TS(1:(7*TS_length));
received_TS_Schmidl=received_TS((TS_length+1):(2*TS_length));
received_TS_Minn=received_TS((3*TS_length+1):(4*TS_length));
received_TS_Park=received_TS((5*TS_length+1):(6*TS_length));
subplot(311)
plot(0:Nt_carr-1+CP_length,received_TS_Schmidl);%接收端Schmidl训练序列波形
axis([0,Nt_carr-1+CP_length, -1.5, 1.5]);
grid on;
ylabel('Amplitude');
xlabel('Time');
title('接收端Schmidl训练序列波形');

subplot(312)
plot(0:Nt_carr-1+CP_length,received_TS_Minn);%接收端Minn训练序列波形
axis([0,Nt_carr-1+CP_length, -1.5, 1.5]);
grid on;
ylabel('Amplitude');
xlabel('Time');
title('接收端Minn训练序列波形');

subplot(313)
plot(0:Nt_carr-1+CP_length,received_TS_Park);%接收端Park训练序列波形
axis([0,Nt_carr-1+CP_length, -1.5, 1.5]);
grid on;
ylabel('Amplitude');
xlabel('Time');
title('接收端Park训练序列波形');

%=====================Schmidl方法进行符号定时====================================
Syn_length=3*TS_length/2;
p_Schmidl=zeros(1,Syn_length);
R_Schmidl=zeros(1,Syn_length);
M_Schmidl=zeros(1,Syn_length);
Syn_Schmidl=received_TS(1:(3*TS_length));
for m=1+CP_length:Syn_length+CP_length
    for k=0:Nt_carr/2-1
        p_Schmidl(m-CP_length)=p_Schmidl(m-CP_length)+conj(Syn_Schmidl(m+k))*Syn_Schmidl(m+k+Nt_carr/2);
        R_Schmidl(m-CP_length)=R_Schmidl(m-CP_length)+abs(Syn_Schmidl(m+k+Nt_carr/2))^2;
    end
    M_Schmidl(m-CP_length)=abs(p_Schmidl(m-CP_length))^2/R_Schmidl(m-CP_length)^2;
end
[maximum_Schmidl,max_index_Schmidl]=max(M_Schmidl);

%=====================Minn方法进行符号定时====================================
p_Minn=zeros(1,Syn_length);
R_Minn=zeros(1,Syn_length);
M_Minn=zeros(1,Syn_length);
Syn_Minn=received_TS((2*TS_length+1):(5*TS_length));
for m=CP_length+1:Syn_length+CP_length
    for l=0:1
        for k=0:Nt_carr/4-1
            p_Minn(m-CP_length)=p_Minn(m-CP_length)+conj(Syn_Minn(m+k+l*Nt_carr/2))*Syn_Minn(m+k+l*Nt_carr/2+Nt_carr/4);
            R_Minn(m-CP_length)=R_Minn(m-CP_length)+abs(Syn_Minn(m+k+l*Nt_carr/2+Nt_carr/4))^2;
        end
    end
    M_Minn(m-CP_length)=abs(p_Minn(m-CP_length))^2/R_Minn(m-CP_length)^2;
end
[maximum_Minn, max_index_Minn]=max(M_Minn);

%=====================Park方法进行符号定时====================================
p_Park=zeros(1,Syn_length);
R_Park=zeros(1,Syn_length);
M_Park=zeros(1,Syn_length);
Syn_Park=received_TS((4*TS_length+1):(7*TS_length));
for m=CP_length+1:3*TS_length/2+CP_length
    for k=0:Nt_carr/2-1
        p_Park(m-CP_length)=p_Park(m-CP_length)+Syn_Park(m+k)*Syn_Park(m-k+Nt_carr-1);
        R_Park(m-CP_length)=R_Park(m-CP_length)+abs(Syn_Park(m+k))^2;
    end
    M_Park(m-CP_length)=abs(p_Park(m-CP_length))^2/R_Park(m-CP_length)^2;
end
[maximum_Park, max_index_Park]=max(M_Park);

figure;
plot(1:Syn_length,M_Schmidl,'r')
hold on
plot(1:Syn_length,M_Minn,'g')
hold on
plot(1:Syn_length,M_Park,'b')
grid on; 
% axis([0,400,0,1.1]); 
title('三种方法得到的定时量度曲线'); 
legend('Schmidl算法','Minn算法','Park算法')
xlabel('Time(sample)'); 
ylabel('Timing Metric');

%===========================串并转换========================================
received_time_wave_matrix_add_CP=reshape(received_time_wave_sequence_add_TS((7*TS_length+1):end)',Nt_carr+CP_length,Sig_per_carr)';%去除训练序列
figure;
subplot(311)
plot(0:Nt_carr-1+CP_length,received_time_wave_matrix_add_CP(1,:));%接收端第一个OFDM符号的时域波形（含CP）
axis([0,Nt_carr-1+CP_length, -0.6, 0.6]);
grid on;
ylabel('Amplitude');
xlabel('Time');
title('接收端第一个OFDM符号的时域波形（含CP）');

%=====================去除循环前缀CP====================================
received_time_wave_matrix=received_time_wave_matrix_add_CP(:,CP_length+1:CP_length+Nt_carr);%去CP
subplot(312)
plot(0:Nt_carr-1,real(received_time_wave_matrix(1,:)));%接收端第一个OFDM符号的时域波形（含CP）
axis([0,Nt_carr-1+CP_length, -0.6, 0.6]);
grid on;
ylabel('Amplitude');
xlabel('Time');
title('接收端去除CP后的第一个OFDM符号的时域波形');

%=====================信道估计与均衡====================================
% TS_FFT=fft(TS_Schmidl((CP_length+1):TS_length));
% received_TS_FFT=fft(received_TS_Schmidl((CP_length+1):TS_length));
% H=received_TS_FFT./TS_FFT;

H=received_TS_Schmidl((CP_length+1):TS_length)./TS_Schmidl((CP_length+1):TS_length);

% H=fft(received_time_wave_matrix(1,:))./fft(time_wave_matrix(1,:));

received_time_wave_matrix_FFT=fft(received_time_wave_matrix,Nt_carr,2);
received_time_wave_matrix_equilibrium_FFT=received_time_wave_matrix_FFT./H;

received_time_wave_matrix_equilibrium=ifft(received_time_wave_matrix_equilibrium_FFT,Nt_carr,2);
subplot(313)
plot(0:Nt_carr-1,real(received_time_wave_matrix_equilibrium(1,:)));%接收端第一个OFDM符号的时域波形（含CP）
axis([0,Nt_carr-1+CP_length, -0.6, 0.6]);
grid on;
ylabel('Amplitude');
xlabel('Time');
title('接收端均衡后第一个OFDM符号的时域波形');

%=================================FFT======================================
received_equilibrium_Hermite=fft(received_time_wave_matrix_equilibrium,Nt_carr,2);

%===========================解埃尔米特映射==================================
received_equilibrium=received_equilibrium_Hermite(:,carriers);

%===========================降低PAPR后处理==================================
received_equilibrium(:,3:Np_carr)=received_equilibrium(:,3:Np_carr)/b(1);
received_equilibrium(:,1)=received_equilibrium(:,1)/b(ceil(index/2));
received_equilibrium(:,2)=received_equilibrium(:,2)/b(index-2*(ceil(index/2)-1));

%=============================并串转换=====================================
received_complex_carrier_matrix1=reshape(received_equilibrium',Np_carr*Sig_per_carr,1)';

figure;
plot(received_complex_carrier_matrix1,'*r');%接收端星座图
title('接收端(16QAM解调制前)星座图');
axis([-5,5,-5,5]);
grid on   %显示网格线
axis square

%===========================16QAM解调==================================


   demodu_baseband_out=qamdemod(received_complex_carrier_matrix1,M,'Gray','UnitAveragePower',true);
    rxBits = de2bi(demodu_baseband_out.','left-msb');
   rxBits2 = reshape(rxBits.',1,baseband_out_length);

[~,ber]=symerr(rxBits2,baseband_out);
ber_carriers=zeros(1,Np_carr);
for j=1:Np_carr
    for i=1:Sig_per_carr
        for k=1:K
            if rxBits2((i-1)*Np_carr*K+(j-1)*K+k)~=baseband_out((i-1)*Np_carr*K+(j-1)*K+k)
                ber_carriers(j)=ber_carriers(j)+1;
            end
        end
    end
end
figure;
plot(1:Np_carr,ber_carriers,'--r*');
title('各子载波误码个数')
ylabel('误码个数');
xlabel('子载波编号');
xlim([1 Np_carr]);
grid on;

部分实验图如下： 

 

 

 

 

标签：received,OFDM,carr,TS,length,MATLAB,16QAM,CP,Nt
来源： https://blog.csdn.net/naturly/article/details/120473814

本站声明： 1. iCode9 技术分享网（下文简称本站）提供的所有内容，仅供技术学习、探讨和分享；
2. 关于本站的所有留言、评论、转载及引用，纯属内容发起人的个人观点，与本站观点和立场无关；
3. 关于本站的所有言论和文字，纯属内容发起人的个人观点，与本站观点和立场无关；
4. 本站文章均是网友提供，不完全保证技术分享内容的完整性、准确性、时效性、风险性和版权归属；如您发现该文章侵犯了您的权益，可联系我们第一时间进行删除；
5. 本站为非盈利性的个人网站，所有内容不会用来进行牟利，也不会利用任何形式的广告来间接获益，纯粹是为了广大技术爱好者提供技术内容和技术思想的分享性交流网站。

关于我们 | 联系我们 | 留言反馈

专注分享技术，共同学习，共同进步。侵权联系[81616952@qq.com]

Copyright (C)ICode9.com, All Rights Reserved.

ICode9版权所有