%OFDM Channel Estimation
%IFFT_bin_length: IFFT和FFT的点数
%carrier_count: 子载波个数
%bits_per_symbol: 每符号上的比特数
%symbols_per_carrier: 每载波的OFDM符号数
%X:欲发送的二进制比特流
clear all;
clc;
IFFT_bin_length=512; %ifft长度
carrier_count=100; %子载波数
bits_per_symbol=2; %每符号比特数,采用4进制调制
symbols_per_carrier=12; %一桢符号数
LI=7; %导频之间的间隔
Np=ceil(carrier_count/LI)+1; %导频数加1的原因:使最后一列也是导频
N_number=carrier_count*symbols_per_carrier*bits_per_symbol; %一祯比特数
carriers=1:carrier_count+Np; %子载波加导频
GI=8; % guard interval length
N_snr=40; % 每比特信噪比
snr=8; % 信噪比间隔
%------------------------------------------------------------
% vector initialization
X=zeros(1,N_number); %2400个bit
X1=[];
X2=[];
X3=[];
X4=[];
X5=[];
X6=[];
X7=[];
Y1=[];
Y2=[];
Y3=[];
Y4=[];
Y5=[];
Y6=[];
Y7=[];
XX=zeros(1,N_number); %2400
dif_bit=zeros(1,N_number); %2400
dif_bit1=zeros(1,N_number); %2400
dif_bit2=zeros(1,N_number); %2400
dif_bit3=zeros(1,N_number); %2400
X=randint(1,N_number); %产生二进制随即序列(非0即1)2400,产生2400个二进制随即序列
%--------------------------------------------------------
%QPSK调制:(1 1)->pi/4;(0 1)->3*pi/4;(0 0)->-3*pi/4;(1,0)->-pi/4;
s=(X.*2-1)/sqrt(2);
sreal=s(1:2:N_number);
simage=s(2:2:N_number);
X1=sreal+j.*simage; %已调信号bit流,共有1200个已调制随即序列的输出。
%---------------------------------------------------------
%产生随机导频信号
%--------------------------------------------------------
train_sym=randint(1,2*symbols_per_carrier); %24个随机数
t=(train_sym.*2-1)/sqrt(2);
treal=t(1:2:2*symbols_per_carrier);
timage=t(2:2:2*symbols_per_carrier);
training_symbols1=treal+j.*timage; % 0.7071 - 0.7071i -0.7071 - 0.7071i -0.7071 - 0.7071i。。。1*12
training_symbols2=training_symbols1.';
training_symbols=repmat(training_symbols2,1,Np); %12*16 复制第一列变成16列
%disp(training_symbols)
pilot=1:LI+1:carrier_count+Np; %导频插入位置序号1 9 17 25 33 41 49 57 65 73 81 89 97 105 113
if length(pilot)~=Np
pilot=[pilot,carrier_count+Np]; %最后一列变成导频1 9 17 25 33 41 49 57 65 73 81 89 97 105 113 116
% pilot=[1 9 17 25 33 41 49 57 65 73 81 89 97 105 113 116]
end
%--------------------------------------------------------
%串并转换
X2=reshape(X1,carrier_count,symbols_per_carrier).'; %12个复信号符号,100列载波
%---------------------------------------------------------
%插入导频
signal=1:carrier_count+Np; % 1:116
signal(pilot)=[]; % 从1:116中挖去了pilot部分
X3(:,pilot)=training_symbols;
% 将training_symbols的第1列赋给X3矩阵的第1列,把training_symbols的第2列赋给X3矩阵的第9列,依此类推
X3(:,signal)=X2; % 再放入12*100,100列子载波,共12*116
% 所得到的X3就是插入了导频之后的并且经过串并转换了的已调制序列。
% X3=cat(1,training_symbols,X2);
FFT_modulation=zeros(symbols_per_carrier,128); %12*128的全0矩阵
FFT_modulation(:,carriers)=X3; %12*128列的矩阵,其中的前116列安放了之前的12*116个数据
FFT_X3=fft(FFT_modulation,128,2);
for i=1:12
inputSamples_ifdma(i,1:4:IFFT_bin_length) = FFT_X3(i,:); % 子载波映射
end
X4=ifft(inputSamples_ifdma,IFFT_bin_length,2); %每个符号512点ifft,最后形成了X4是12*512的经过IFFT之后的矩阵
%加保护间隔(循环前缀)
for k=1:symbols_per_carrier; % 一共有12行,所以要做12次循环
for i=1:IFFT_bin_length;
X6(k,i+GI)=X4(k,i);
end
for i=1:GI;
X6(k,i)=X4(k,i+IFFT_bin_length-GI);
end
end
%---------------------------------------------------------
%并串转换
X7=reshape(X6.',1,symbols_per_carrier*(IFFT_bin_length+GI));
%12*520先转置,再变成1*6240复信号流
%---------------------------------------------------------
%信道模型:带多普勒频移的瑞利衰落信道
fd=300; %多普勒频移
r=6; %多径数
a=[0.123 0.3 0.4 0.5 0.7 0.8]; %多径的幅度
d=[2 3 4 5 9 13]; %各径的延迟
T=1; %系统采样周期
th=[90 0 72 144 216 288]*pi./180; %相移
h=zeros(1,carrier_count); %1*100
hh=[];
for k=1:r %1:6
%deta=[zeros(1,d(k)-1),1,zeros(1,carrier_count-d(k))];
h1=a(k)*exp(j*((2*pi*T*fd*d(k)/carrier_count)));
%h1=a(k)*exp(j*((2*pi*T*fd*d(k)/carrier_count)));
hh=[hh,h1];
end
h(d+1)=hh; % 3 4 5 6 10 14 处有多径效应
%noise=randn(1,length(X7))+j.*randn(1,length(X7));
%--------------------------------------------------------
channel1=zeros(size(X7)); % 1*1632
channel1(1+d(1):length(X7))=hh(1)*X7(1:length(X7)-d(1));
channel2=zeros(size(X7));
channel2(1+d(2):length(X7))=hh(2)*X7(1:length(X7)-d(2));
channel3=zeros(size(X7));
channel3(1+d(3):length(X7))=hh(3)*X7(1:length(X7)-d(3));
channel4=zeros(size(X7));
channel4(1+d(4):length(X7))=hh(4)*X7(1:length(X7)-d(4));
channel5=zeros(size(X7));
channel5(1+d(5):length(X7))=hh(5)*X7(1:length(X7)-d(5));
channel6=zeros(size(X7));
channel6(1+d(6):length(X7))=hh(6)*X7(1:length(X7)-d(6));
%---------------------------------------------------------------
Tx_data=X7+channel1+channel2+channel3+channel4; %4径干扰后的数据流1*1632
%---------------------------------------------------------------
%---------------------------------------------------------------
%----------------------------------------------------------------
%加高斯白噪声
Error_ber=[]; %误比特率
Error_ber1=[];
Error_ber2=[]; %误比特率
Error_ber3=[];
%Error_ser=[]; %误符号率
for snr_db=0:snr:N_snr %0:8:40共6个数
code_power=0;
code_power=[norm(Tx_data)]^2/(length(Tx_data)); %信号的符号功率
%bit_power=var(Tx_data);
bit_power=code_power/bits_per_symbol; %比特功率
noise_power=10*log10((bit_power/(10^(snr_db/10)))); %噪声功率
noise=wgn(1,length(Tx_data),noise_power,'complex'); %产生GAUSS白噪声信号
Y7=Tx_data+noise;
%-------------------------------------------------------
%串并变换
Y6=reshape(Y7,IFFT_bin_length+GI,symbols_per_carrier).';%先变成520*12,再转置成12*520,恢复成行为符号,列为载波
%去保护间隔
for k=1:symbols_per_carrier;
for i=1:IFFT_bin_length;
Y5(k,i)=Y6(k,i+GI);
end
end
Y4 = fft(Y5,IFFT_bin_length,2); % 每行的符号进行512点fft 12*512
for j=1:12
Y4_1(j,:) = Y4(j,1:4:IFFT_bin_length); % 恢复子载波映射 Y3是12*128的矩阵
end
Y3_1 = ifft(Y4_1,128,2); % 对12*128的矩阵对每行的IFFT
Y3 = Y3_1(:,1:116);
%-------------------------------------------------------------
%LS信道估计
%H=[];
Y2=Y3(:,signal); % 将全部的数据列取出来,一共是12行100列
Rx_training_symbols=Y3(:,pilot);
length(Rx_training_symbols)
%上述两条语句就是将信息序列的位置和训练导频序列的位置进行了分离
Rx_training_symbols0=reshape(Rx_training_symbols,symbols_per_carrier*Np,1);
training_symbol0=reshape(training_symbols,1,symbols_per_carrier*Np);
training_symbol1=diag(training_symbol0);
%disp(training_symbols)