echo off;clear all;close all;clc;
fprintf( 'OFDM仿真\n') ;
tic
% --------------------------------------------- %
% Parameter Declaration %
% --------------------------------------------- %
% Initialize the parameters
NumLoop = 1000;
NumSubc = 128;
NumCP = 8;
SyncDelay = 0;
%-----------------------------------------------
% 子载波数 128
% 位数/ 符号 2
% 帧个数 1000
% 训练符号数 0
% 循环前缀长度 8 (1/16)*T
% 调制方式 BPSK
% 多径信道数 5
% IFFT Size 128
% 信道最大时延 2
% --------------------------------------------- %
% BPSK Modulation %
% --------------------------------------------- %
% Generate the random binary stream for transmit test
BitsTx = floor(rand(1,NumLoop*NumSubc)*2);
BPSKTable = [-1 1];
SymBPSK = BPSKTable(BitsTx+1); % Modulate (Generates BPSK symbols)
% --------------------------------------------- %
% IFFT %
% --------------------------------------------- %
SymIFFT = zeros(NumSubc,NumLoop);
SymIFFTtmp = reshape(SymBPSK,NumSubc,NumLoop);
SymIFFT = ifft(SymIFFTtmp,NumSubc,1);
% --------------------------------------------- %
% Add cyclic prefix %
% --------------------------------------------- %
NumAddPrefix = NumSubc + NumCP;
SymCP = zeros(NumAddPrefix,NumLoop);
RowPrefix = (NumSubc - NumCP + 1):NumSubc;
SymCP = [SymIFFT(RowPrefix,:);SymIFFT];
% --------------------------------------------- %
% Go through the channel %
% --------------------------------------------- %
SymCh = zeros(1,NumAddPrefix*NumLoop);
SymChtmp = SymCP(:).'; % 进行这个转置操作之后就成了一个矢量
% 相当于把矩阵的列向量依次排列 改变为一个行向量
Ch = [0.227 0.46 0.688 0.46 0.227];
SymChtmptmp = filter(Ch,1,SymChtmp);
%--------------------------------------------------------------------------
% 函数说明：
% Firlter data with an infinite impulse response (IIR) or finite impulse response
% (FIR) filter
% y = filter(b,a,X) filters the data in vector X with the filter described by
% numerator coefficient vector b and denominator coefficient vector a. If a(1) is
% not equal to 1, filter normalizes the filter coefficients by a(1). If a(1) equals
% 0, filter returns an error.
%--------------------------------------------------------------------------
% If X is a matrix, filter operates on the columns of X. If X is a multidimensional
% array, filter operates on the first nonsingleton dimension.
%-------------------------------------------------------------------------
BerSnrTable = zeros(30,5);
for snr=0:49;
BerSnrTable(snr+1,1) = snr;
BerSnrTable(snr+1,2) = 10^(BerSnrTable(snr+1,1)/10); %将信噪比由db变换到一般数据比
BerSnrTable(snr+1,3) = 1/sqrt(2*BerSnrTable(snr+1,2)); %信道的噪声标准差
SymCh = AWGN(SymChtmptmp,snr,'measured'); % Add the noise
%--------------------------------------------------------------------------
% 函数说明：
% AWGN Add white Gaussian noise to a signal.
% Y = AWGN(X,SNR) adds white Gaussian noise to X. The SNR is in dB.
% The power of X is assumed to be 0 dBW. If X is complex, then
% AWGN adds complex noise.
% ------------------------------------------------------------------------
% Y = AWGN(X,SNR,SIGPOWER) when SIGPOWER is numeric, it represents
% the signal power in dBW. When SIGPOWER is 'measured', AWGN measures
% the signal power before adding noise.
% -------------------------------------------------------------------------
% Y = AWGN(X,SNR,SIGPOWER,STATE) resets the state of RANDN to STATE.
%
% Y = AWGN(..., POWERTYPE) specifies the units of SNR and SIGPOWER.
% POWERTYPE can be 'db' or 'linear'. If POWERTYPE is 'db', then SNR
% is measured in dB and SIGPOWER is measured in dBW. If POWERTYPE is
% 'linear', then SNR is measured as a ratio and SIGPOWER is measured
% in Watts.
%
% Example: To specify the power of X to be 0 dBW and add noise to produce
% an SNR of 10dB, use:
% X = sqrt(2)*sin(0:pi/8:6*pi);
% Y = AWGN(X,10,0);
%
% Example: To specify the power of X to be 0 dBW, set RANDN to the 1234th
% state and add noise to produce an SNR of 10dB, use:
% X = sqrt(2)*sin(0:pi/8:6*pi);
% Y = AWGN(X,10,0,1234);
%
% Example: To specify the power of X to be 3 Watts and add noise to
% produce a linear SNR of 4, use:
% X = sqrt(2)*sin(0:pi/8:6*pi);
% Y = AWGN(X,4,3,'linear');
%
% Example: To cause AWGN to measure the power of X, set RANDN to the
% 1234th state and add noise to produce a linear SNR of 4, use:
% X = sqrt(2)*sin(0:pi/8:6*pi);
% Y = AWGN(X,4,'measured',1234,'linear');
% --------------------------------------------- %
% Remove Guard Intervals %
% --------------------------------------------- %
SymDeCP = zeros(NumSubc,NumLoop);
SymDeCPtmp = reshape(SymCh,NumSubc + NumCP,NumLoop);
SymDeCP = SymDeCPtmp((NumCP+1+SyncDelay):NumAddPrefix+SyncDelay,:);
% --------------------------------------------- %
% FFT %
% --------------------------------------------- %
SymFFTtmp = fft(SymDeCP,NumSubc,1);
% --------------------------------------------- %
% Equalize %
% --------------------------------------------- %
H = fft(Ch,NumSubc);
for i=1:NumSubc
W(i) = conj(H(i))/(BerSnrTable(snr+1,3)^2+abs(H(i))^2);
end
E_W=zeros(NumSubc,NumLoop);
for n1=1:NumLoop
E_W(:,n1)=W;
end
SymFFT = SymFFTtmp.*E_W;
% --------------------------------------------- %
% Make Decision(Include DeBPSK) %
% --------------------------------------------- %
BitsRx = zeros(1,NumSubc*NumLoop);
BitsRxtmp = SymFFT(:).';
for m = 1:NumLoop*NumSubc
Real =real(BitsRxtmp(1,m));
if Real>0
BitsRx(1,m)=1;
else
BitsRx(1,m)=0;
end
end
% --------------------------------------------- %
% Compute Error Number/Rate %
% --------------------------------------------- %
[Num,Ber] = symerr(BitsTx,BitsRx)
BerSnrTable(snr+1,4) = Num ;
BerSnrTable(snr+1,5) = Ber ;
end
figure(1);
subplot(2,1,1);
semilogy(BerSnrTable(:,1),BerSnrTable(:,4),'o-');
subplot(2,1,2);
semilogy(BerSnrTable(:,1),BerSnrTable(:,5),'o-');
% ------------------------------------------------------------------------
time_of_sim = toc
echo on;
% --------------------------------------------- %
% The END %
% --------------------------------------------- %
