cp_ofdm_code 做ofdm的matlab仿真代码 修改一下文件夹路径即可运行 可以作为xuexizili

SHANHONGBO

% clear all;
%---------------------------------------仿真参数的设定--------------------------------------------%
N = 48;                                    %有效子载波个数
N_number = 64;                             %子载波个数
N_symbol = 20;                             %symbol个数
CP_Length = N_number/4;                    %CP(Cyclic Prefix)的长度
Symbol_Length = CP_Length + N_number;      %每个符号的长度
MQAM = 4;                                  %4QAM调制
pilot_symbol_number = 2;                   %导频符号个数


%---------------------------------------产生数据符号---------------------------------------------%
% dsource_buf = zeros(1,N*N_symbol);
%
% dsource_buf = randn(1,2*N*N_symbol);
% dsource_buf = dsource_buf > 0;                       %转换成0、1随机序列
load dsource_buf.mat

%----------------------------扰码---------------------------------%
scrambler_out=zeros(1,1920);
reg = [1, 0, 1, 1, 1, 0, 1];
for i = 1: 1920                                          
   scrambler_out(1,i) = mod(reg(7) + reg(4) + dsource_buf(1,i) , 2);   % g(x) = x^7 + x^4 + 1
   x_1 = mod(reg(7) + reg(4) , 2);                       % x^1 =  x^7 + x^4     
   reg = [x_1 , reg(1:end-1)];                            % register refresh
end

scrambler_out = [scrambler_out(1:1914),0,0,0,0,0,0];%生成一帧内的随机二进制比特流，此帧为发送端发送的数据

%-----------------------------------------FEC coder----------------------------------------------%
% coder_constraint = 7;
% d_trellis = poly2trellis(7,[133 171]);               %设置卷积码的参数，码率为1/2，约束长度为7的卷积码
% out_fec = convenc(scrambler_out,d_trellis);                %使用系统函数完成卷积码的编码

%----------interleaving----------
%
% k=0:2*N*N_symbol-1;                   %first step interleaving
% s=(2*N*N_symbol/16)*mod(k,16)+floor(k/16);
% out_interleaving(k+1)=out_fec(s+1);
%   
%----------------------------------------OQAM mapping------------------------------------------------%
data = zeros(N_number,N_symbol);

out_map_i = 2*scrambler_out(1:2:2*N*N_symbol-1)-1;        %把奇数位上的bit进行映射，0—>-1,1—>1
out_map_q = 2*scrambler_out(2:2:2*N*N_symbol)-1;          %把偶数位上的bit进行映射,0—>-1,1—>1
out_map = out_map_i + 1i*out_map_q;                 %把奇数位的映射作为实部，把偶数位的映射作为虚部
% for m=1:N_symbol
%     data(1:N,m)=out_map(N*(m-1)+1:N*m);
% end
for m = 1:N_symbol                                  %数据部分分布在48个有效子载波上
    data(2:7 , m) =  out_map(1 , 48*(m - 1) + 1 : 48*(m-1)+6);
    data(9:21 , m) =  out_map(1 , 48*(m - 1) + 7 : 48*(m - 1)+19);
    data(23:27 , m) = out_map(1 , 48*(m - 1) + 20 : 48*(m-1)+24);
    data(39:43 , m) = out_map(1 , 48*(m - 1) + 25 : 48*(m-1)+29);
    data(45:57 , m) = out_map(1 , 48*(m - 1) + 30 : 48*(m-1)+42);
    data(59:64 , m) = out_map(1 , 48*(m - 1) + 43 : 48*m);
    data(8,m) = 1;                                  %在数据中插入4个梳状导频，共52个子载波
    data(22,m) = -1;
    data(44,m) = 1;
    data(58,m) = 1;
end

pilot1 = [0,1, -1, -1, 1, 1, -1, 1, -1, 1, -1, -1, -1, -1, -1, 1, 1, -1, -1, 1, -1, 1, -1, 1, 1, 1, 1,0,0,0,0,0,0,0,0,0,0,0,1, 1, -1, -1, 1, 1, -1, 1, -1, 1, 1, 1, 1, 1, 1, -1, -1, 1, 1, -1, 1, -1, 1, 1, 1, 1].';
pilot = [pilot1,pilot1];
data_final = [pilot,data];
%--------------------------------------------------------IFFT的实现方式-------------------------------------------------------------%
signal_ifft = zeros(1,(N_symbol+pilot_symbol_number)*Symbol_Length);     %信号长度为（CP长度+IFFT长度）乘以符号个数

for m = 1N_symbol+pilot_symbol_number)
    x = ifft(data_final(1:N_number,m));           %IFFT
    x = x.';                                %将列向量转换成行向量，并串转换
    signal_ifft(1,(m-1)*Symbol_Length+1:m*Symbol_Length) = [x(N_number-CP_Length+1:N_number) , x];    %加入CP
end
signal_transmit = signal_ifft;

imag_max = max(imag(signal_transmit));
real_max = max(real(signal_transmit));
IQ_max = max(imag_max,real_max);

signal_transmit = signal_transmit ./ IQ_max .* 2^15;

Q_data=imag(signal_transmit);
I_data=real(signal_transmit);

Q_data = [zeros(1,400) ,repmat(Q_data,1,10)];
I_data = [zeros(1,400) ,repmat(I_data,1,10)];


fp=fopen('E:\shb\Q_data.txt','w');
fprintf(fp,'%d\n',Q_data);
fclose(fp);

fp=fopen('E:\shb\I_data.txt','w');
fprintf(fp,'%d\n',I_data);
fclose(fp);

% ------------------------接收端---------------------------%
%  signal_received = signal_transmit;signal_received
Rxdata= GetRxData('E:\shb\RX2.bin',10000);

plot(abs(Rxdata));

signal_received = Rxdata(385+1 : 385 + 1760);

% fs=2*e^6;
% dt=1/fs;
% [c,lags]=xcorr(signal_transmit,signal_received);
% [~,time]=max(abs(c));
% delay=lags(time);
% delaytime=delay./fs;
% figure
% plot(lags,c);
%--------------------------------------------------------FFT的方式解调-------------------------------------------------------------%
rece_data = zeros(N_number,N_symbol+pilot_symbol_number);

for m = 1:N_symbol+pilot_symbol_number
    rec_X = signal_received(1,(m-1)*Symbol_Length+1:m*Symbol_Length);           %取出每帧符号
    rece_data(1:N_number,m) = fft(rec_X(CP_Length+1:Symbol_Length));         %去掉CP，并且做FFT变换
end


%----------------------------------------信道估计-------------------------------------------------%
H_channel1 = zeros(N_number,1);
H_channel1(2:27) = rece_data(2:27,1)./data_final(2:27,1);
H_channel1(39:64) = rece_data(39:64,1)./data_final(39:64,1);

H_channel2 = zeros(N_number,1);
H_channel2(2:27) = rece_data(2:27,2)./data_final(2:27,2);
H_channel2(39:64) = rece_data(39:64,2)./data_final(39:64,2);

H_channel = (H_channel1+H_channel2)/2;

Hmn_equalizer = zeros(N_number,1);
Hmn_equalizer(2:27) = 1 ./  H_channel(2:27);
Hmn_equalizer(39:64) = 1 ./ H_channel(39:64);

%%----------------------------------用CSI做均衡------------------------------%%
signal_final = zeros(N_number,N_symbol);
for i = 1 : N_symbol
    signal_final(:,i) =  rece_data(:,i+pilot_symbol_number).* Hmn_equalizer;
     %   signal_final(:,i) =  rec_data(:,i+pilot_symbol_number);
end

% %--------------------------------计算信噪比---------------------------------%
noise= zeros(N,N_symbol);
noise(1:6,:)=data(2:7,:)-signal_final(2:7,:);
noise(7:19,:)=data(9:21,:)-signal_final(9:21,:);
noise(20:24,:)=data(23:27,:)-signal_final(23:27,:);
noise(25:29,:)=data(39:43,:)-signal_final(39:43,:);
noise(30:42,:)=data(45:57,:)-signal_final(45:57,:);
noise(43:48,:)=data(59:64,:)-signal_final(59:64,:);
noise_power=noise.*conj(noise);
signal_n=sum(noise_power,2);
signal_m=sum(signal_n,1);
SNR = 10*log10(2*N_symbol*N/signal_m);


%------------------将中间插零的64行信号 变为48行数据-----------------------------

signal_final_data=zeros(N,N_symbol);
signal_final_data(1:6,:)=signal_final(2:7,:);
signal_final_data(7:19,:)=signal_final(9:21,:);
signal_final_data(20:24,:)=signal_final(23:27,:);
signal_final_data(25:29,:)=signal_final(39:43,:);
signal_final_data(30:42,:)=signal_final(45:57,:);
signal_final_data(43:48,:)=signal_final(59:64,:);

%---------------------------------画星座图----------------------------------%
signal_final_data_star_map = zeros(1,N*N_symbol);
for m=1:N_symbol
    signal_final_data_star_map(N*(m-1)+1:N*m)=signal_final_data(1:N,m);
end
scatterplot(signal_final_data_star_map(1,:));



%----------------------------------------demapping------------------------------------------------%
out_demapping = zeros(1,2*N*N_symbol);

for m = 1:N_symbol
    out_demapping(1,2*N*(m-1)+1:2:2*N*m) = real(signal_final_data(1:N,m));
    out_demapping(1,2*N*(m-1)+2:2:2*N*m) = imag(signal_final_data(1:N,m));
end
out_demapping = out_demapping > 0;

%----------deinterleaving解交织-----------------%
% out_deinterleaving = zeros(1,N*N_symbol*2);
% j=0:N*N_symbol*2-1;
% k=(N*N_symbol*2/16)*j-(N*N_symbol*2-1)*floor(j/16);
% out_deinterleaving(k+1) = out_demapping(j+1);

%-----------------------------------------FEC decoder------------------------------------------------%
% out_decoder = zeros(1,2*N*N_symbol);
% out_decoder = vitdec(out_deinterleaving,d_trellis,coder_constraint*4,'term','hard');


%----------解扰码------------------------------------------------------%
descrambler_out=zeros(1,1920);
reg = [1, 0, 1, 1, 1, 0, 1];
for i = 1: 1920                                          
   descrambler_out(1,i) = mod(reg(7) + reg(4) + out_demapping(1,i) , 2);   % g(x) = x^7 + x^4 + 1
   x_1 = mod(reg(7) + reg(4) , 2);                       % x^1 =  x^7 + x^4     
   reg = [x_1 , reg(1:end-1)];                            % register refresh
end

BER = sum(abs(descrambler_out(1:end-6)-dsource_buf(1:end-6))) / length(dsource_buf(1:end-6));